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Figure  21.  Responses  as  functions  of  frequency  and  rms  free  surface  27 

fluctuations.  Top  row:  the  near  source  responses  of  direct  Rg ,  scattered  Rg  and 
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Figure  22.  Similar  to  Figure  4  except  each  panel  is  for  a  different  source  depth.  29 

Figure  23.  Similar  to  Figure  5,  except  responses  are  functions  of  frequencies  and 
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1.  SUMMARY 


The  complex  excitation  and  energy  partitioning  mechanisms  yielding  regional 
phases  are  difficult  to  empirically  separate  by  data  analysis.  Thus,  numerical  modeling 
approaches  are  valuable  for  investigating  excitation  and  propagation  of  regional  seismic 
phases.  We  use  accurate  full-wave  simulations  (2D  and  3D  finite-difference  method,  and 
2D  boundary  element  method)  to  calculate  seismic  wave  excitation  and  propagation  in 
the  near-source  region.  An  embedded  array  slowness  analysis  is  used  for  quantifying  how 
wavefield  energy  will  be  partitioned  into  the  long-range  propagation  regime.  Due  to  its 
high  efficiency,  the  method  can  simulate  near-source  processes  using  very  fine  structures. 
A  large  number  of  source  and  model  parameters  can  be  examined  for  broad  frequency 
ranges.  We  use  this  method  to  investigate  the  effect  of  volumetric  and  topographic 
scattering  on  the  near-source  energy  partitioning  for  an  explosion  source.  Different 
random  velocity  fluctuations  and  topographic  models,  variable  source  depth,  and 
different  Q  models  are  investigated  using  numerical  simulations.  The  responses  of 
different  phases  as  functions  of  frequency  and  source/model  parameters  are  calculated 
and  their  energy  budget  evaluated. 

The  results  reveal  that  the  depth  of  the  source  and  the  depth  of  the  scattering 
process  have  strong  effects  on  P-to-S  conversion  and  partitioning  of  energy  into  trapped 
or  leaking  signals.  The  modeling  shows  that  S*-to-Lg  excitation  is  generally  stronger  for 
low-frequencies  and  shallow  source  depths  while  P-pS-to-Lg  scattering  is  stronger  for 
high-frequencies.  The  free  surface  scattering  has  strong  effect  on  near  source  energy 
partitioning.  The  scattering  process  can  excite  the  Rg-wave  for  a  moderately  rugged 
topography,  but  prevents  the  formation  and  propagation  of  short  period  Rg- waves  when 
the  surface  becomes  too  rugged.  For  models  with  a  high  velocity  shallow  crust,  the  free 
surface  scattering  provides  an  important  mechanism  that  transfers  energy  for  an 
explosion  source  into  the  Lg-wave  in  the  near  source  region.  At  lower  frequencies  and  for 
a  moderately  rugged  free  surface,  the  Rg-to-Lg  transfer  is  relatively  efficient.  At  higher 
frequencies  and  for  a  very  rugged  free  surface,  the  body-to-Lg  transfer  may  dominate  the 
process.  The  correlation  length  of  the  random  free  surface  fluctuation  imparts  specific 
frequency  dependence  to  the  transfer  function,  with  maximum  coupling  near  ka=l. 
Intrinsic  attenuation  within  the  uppermost  crust  has  a  strong  effect  on  the  energy  transfer 
through  surface  scattering,  with  high  frequency  content  losing  energy  faster  than  the 
lower  frequency  waves. 

We  expand  our  simulation  and  analysis  to  fully  3D  models  with  both  volumetric 
heterogeneity  and  rough  free  surface.  Investigating  the  near  source  process  in  3D  is  more 
realistic  because  it  gives  correct  scattering  pattern  and  geometric  spreading.  In  addition,  a 
general  3D  model  provides  the  coupling  between  the  P-SV  and  SH  waves.  Velocity 
models  with  different  volume  heterogeneities  and  surface  fluctuations  are  tested  and 
compared.  For  an  isotropic  explosion  source  in  a  typical  layered  model  with  high  velocity 
top  layer,  only  weak  Lg  wave  (trapped  guided  wave)  can  be  generated  from  the  S*  phase. 
There  is  no  energy  in  the  tangential  component.  When  adding  volume  heterogeneities  or 
surface  fluctuations  in  the  model,  the  energy  transfers  from  leaky  mode  to  trapped  guided 
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waves  as  well  as  from  P-SV  mode  to  the  SH  mode.  For  the  surfaee  seattering,  an 
exponential  random  surface  fluctuation  provide  stronger  coupling  than  a  Gaussian 
random  surface  model.  Generally,  a  model  with  more  small  scale  topographic 
fluctuations  tends  to  generate  more  trapped  energy. 


2.  INTRODUCTION 

With  the  current  emphasis  on  global  monitoring  for  low-yield  nuclear  tests, 
regional  seismic  phases  such  as  Lg  have  become  very  important  for  magnitude  and  yield 
estimation  of  underground  nuclear  tests,  (e.g.,  Nuttli,  1986;  Xie,  et  ah,  1996;  Patton, 
2001).  In  addition,  various  F’/S'-type  amplitude  ratios  for  high  frequency  regional  phases 
(e.g.,  PnISn,  PnlLg,  PgILg,  Pg/Sn)  have  become  important  for  event  discrimination  (e.g., 
Taylor  et  ah,  1989;  Kim  et  ah,  1993,  1997;  Walter  et  ah,  1995;  Fisk  et  ah,  1996;  Taylor, 
1996;  Taylor  and  Hartse,  1997;  Hartse  et  ah,  1997;  Fan  and  Lay,  1998a-c;  Xie,  2002; 
Bottone  et  ah,  2002).  The  applications  of  regional  phases  for  yield  estimation  and  event 
discrimination  are  largely  based  on  empirical  approaches,  and  while  very  promising  in 
many  cases,  there  are  major  questions  about  the  nature  of  excitation  of  ^-wave  dominated 
phases  such  as  Lg.  There  are  similar  questions  regarding  the  relative  excitation  effects  for 
F’/S'-type  ratios  in  regional  phases,  particularly  given  the  huge  scatter  observed  in  both 
earthquake  and  explosion  data  populations. 

There  are  now  many  observational  and  theoretical  studies  addressing  the  regional 
phase  energy  partitioning  issue  in  both  the  near-source  environment  and  the  propagation- 
path  environment.  Along  the  propagation  path,  the  existence  of  small-scale 
heterogeneities  in  the  crust  and  the  associated  seismic  wave  scattering  has  long  been 
addressed  by  seismologists  (Wu  and  Aki,  1988;  Sato  and  Fehler,  1998),  and  its  effect  on 
the  long-range  energy  partitioning  has  been  documented  (Wu,  et  ah,  2000a,  b).  In  the 
source  region,  several  possible  near-source  energy  excitation  mechanisms  have  been 
proposed,  including  P-to-Lg  scattering,  pS-to-Lg  conversion  at  the  free  surface,  Rg-to-Lg 
coupling,  S’^-to-Lg  conversion,  spall  excitation  of  S,  tectonic  release,  explosion 
asymmetry  (CLVD  component)  and  rock-damage  (e.g..  Day  and  Mclaughlin,  1991; 

Gupta  et  ah,  1992,  1997;  Wallace,  1991;  Gutowski,  et  ah,  1984;  Lilwall,  1988;  Xie  and 
Lay,  1994;  Vogljord,  1997;  Johnson  and  Sammis,  2001). 

Due  to  the  complex  excitation  and  energy  partitioning  processes  associated  with 
regional  phases,  it  is  difficult  to  empirically  separate  the  contribution  of  individual  energy 
partitioning  mechanisms  by  analysis  of  data.  Numerical  modeling  approaches  are  thus  of 
great  importance  for  investigating  the  excitation  and  propagation  of  regional  phases. 
Kennett  and  Mykkeltveit  (1984)  and  Kennett  (1989)  used  the  coupled  mode  method  to 
calculate  Lg-wave  propagation  in  crustal  waveguides  with  weak  lateral  heterogeneities. 
Xie  and  Lay  (1994)  investigated  Lg-wave  excitation  using  the  full-wave  finite-difference 
method.  Jih  (1995,  1996)  investigated  Rg-to-Lg  coupling  as  a  possible  Lg  excitation 
mechanism.  With  an  anelastic  FD  method,  Bradley  and  Jones  (1998,  1999)  investigated 
Lg  propagation  and  attenuation  in  Western  China  and  India.  Recently,  using  the  2D  and 
3D  general  Fourier  methods,  Bonner  et  al.  (2003)  investigated  Rg  and  Lg  generation,  and 
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partially  reproduced  the  observed  spectrum  from  the  Depth  of  Burial  Experiment. 

Stevens  et  al.  (2003,  2005),  and  Myers  et  al.  (2005)  investigated  the  physical  basis  of 
explosion  generated  S'-waves  using  a  2D  nonlinear  FD  method,  which  handles 
axisymmetric  near-source  effects  including  spall,  cracking,  and  nonlinear  deformation. 

The  main  disadvantages  of  these  numerical  methods  are  their  low  computation 
efficiency  and  huge  computer  memory  requirement,  especially  when  applied  to 
investigate  broadband  Lg  excitation.  For  the  purpose  of  small  nuclear  test  monitoring, 
the  range  of  interest  for  Fg-wave  simulation  involves  a  broad  frequency  band  (0.2  to  10 
Hz)  and  long  propagation  distances  (up  to  1000  km  or  more).  At  the  same  time,  factors 
that  control  the  source  energy  partitioning  depend  on  the  detailed  source  mechanism  and 
fine  near-source  velocity  structure.  In  addition,  there  are  multiple  mechanisms  that  may 
potentially  contribute  to  the  energy  partitioning  process.  Numerous  parameters  need  to 
be  tested  to  investigate  the  characteristics  of  these  mechanisms,  especially  their 
contributions  to  the  frequency  dependent  features  of  observable  discriminants.  If  random 
heterogeneities  are  to  be  considered,  as  is  likely  to  be  important  for  high  frequency 
signals,  the  results  have  to  be  calculated  statistically  from  simulations  using  a  large 
number  of  realizations.  These  factors  limit  the  approach  of  complete  FD  synthesis  of 
distant  waveforms  for  actual  recording  geometries. 

Although  there  are  continuing  controversies  about  the  dominant  P-io-S  transfer 
mechanisms  affecting  regional  phases,  most  investigators  agree  that  appreciable  energy 
from  explosion  sources  is  converted  to  ^-waves  in  the  near-source  region  (e.g.,  Myers  et 
ah,  2003).  The  physical  processes  by  which  an  explosion  source  generates  regional 
phases  can  be  described  as  energy  partitioning  taking  place  in  the  near-source  region. 

The  partitioned  energy  subsequently  propagates  through  a  long  waveguide,  where 
secondary  energy  partitioning  effects  may  occur,  but  these  are  less  affected  by  the  type  of 
source  involved.  If  the  propagation  effect  is  not  of  primary  interest,  it  is  desirable  to 
avoid  calculating  the  immensely  time-consuming  long-distance  propagation  part  of  the 
problem.  In  order  to  focus  on  the  near-source  energy  partitioning  processes,  we 
developed  a  method  which  uses  accurate  full-wave  simulation  in  the  near-source  region 
which  is  followed  by  a  local  slowness  analysis.  This  method  investigates  energy 
partitioning  right  at  the  source  region,  but  quantifies  how  energy  will  transfer  into  the 
long  range  propagation  path,  which  is  critical  for  comparisons  with  data.  The  localized 
analysis  thereby  isolates  the  physical  processes  controlling  the  energy  partitioning. 

In  the  following  sections,  we  first  present  the  method  with  numerical  calculations 
used  to  demonstrate  its  validity.  Then,  we  use  this  method  to  investigate  several  potential 
Fg-wave  excitation  mechanisms  involving  scattering.  In  section  4,  by  using  2D  finite- 
difference  simulation  in  the  near  source  region,  we  investigate  the  contribution  from  near 
source  volumetric  scattering.  In  section  5,  with  a  2D  boundary  element  method,  we 
investigate  the  contribution  from  the  rough  topography.  In  section  6,  we  use  a  3D  finite- 
difference  method  to  investigate  the  contribution  from  3D  scattering,  including  both 
volumetric  and  surface  effects. 
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3.  METHODOLOGY 


3.1.  The  Near-Source  Strategy 

We  limit  our  detailed  near-souree  simulation  to  a  relatively  small  model  and 
analyze  the  wavefield  within  the  model  to  investigate  the  souree  energy  partitioning  and 
the  exeitation  of  trapped  regional  phases  sueh  as  Lg.  The  numerieal  methods  used  for  this 
purpose  should  be  eapable  of  handling  eomplex  veloeity  models,  e.g.,  volume 
heterogeneities  and  rough  free-surfaee,  and  different  types  of  sourees.  Here  we  use  the 
full-wave  finite-differenee  or  finite-element  methods  for  2D  ealeulation  and  use  a  finite- 
differenee  method  with  rough  topography  eapability  for  3D  simulation.  Many  authors 
(e.g.  Frankel,  1989;  Xie  and  Lay,  1994;  Vogljord,  1997)  have  pointed  out  that  for  S-wave 
energy  to  be  trapped  in  the  waveguide,  reverberating  to  generate  the  Lg-wave,  it  must 
propagate  with  post  eritieal  angle  at  the  Moho  diseontinuity.  However,  in  the  waveguide 
and  espeeially  in  the  near-souree  region,  the  wavefield  is  highly  eomplex.  It  is 
impraetical  to  traee  eaeh  phase  in  the  spaee-time  domain.  An  alternate,  but  equally  valid 
way  of  traeking  the  wave  energy  is  in  the  slowness  domain.  Multiply-reflected  waves 
may  arrive  simultaneously  in  time,  but  in  the  slowness  domain  their  energy  distribution 
gives  clear  information  about  the  wave  intensity,  slowness,  and  propagation  direction. 
Several  methods  can  be  used  to  transfer  space-time  domain  data  into  slowness  (or 
equivalently  wavenumber)  domain  information,  for  example,  FK  analysis  or  slant 
stacking.  Here  we  employ  a  time-domain  local  slant  stacking  method  (Xie  and  Lay,  1994; 
Xie  et  al.  2005a)  to  conduct  slowness  analysis,  working  simultaneously  in  both  the  space 
and  slowness  domains.  The  method  has  two  major  advantages.  First,  it  allows  analysis  of 
the  near-souree  processes  in  multiple  domains,  including  space,  time,  slowness  and 
frequency.  This  allows  isolation  of  different  excitation  and  partitioning  mechanisms 
within  the  complex  near-souree  environment.  Second,  the  embedded  array  slowness 
analysis  method  can  be  applied  at  a  close  range,  well  before  the  Lg-wave  is  actually 
formed.  This  allows  numerous  simulations  for  a  relatively  small  scale  velocity  model 
with  very  fine  near-souree  structures.  Since  this  is  a  computationally  efficient  approach, 
we  can  investigate  energy  partitioning  mechanisms  over  a  broad  frequency  range  and  a 
large  source/model  parameter  space.  For  models  with  random  structural  features,  multiple 
realizations  (hundreds)  are  calculated  and  the  results  are  treated  statistically. 

3.2.  Local  Slowness  Analysis 

Since  the  time-domain  wavefield  does  not  directly  provide  the  propagation 
direction,  we  can  decompose  the  wavefield  into  vector  slowness  domain  using  a  local 
slant  stacking  method  (Xie  and  Lay,  1994;  Xie  et  ah,  2005a) 

u(r,p,l,®)  =  -  r)u[r',l  -  p  •  (r'  -  r),®] ,  (1) 

r' 

where  r  is  the  3D  position  vector,  r'is  the  location  of  the  receiver,  t  is  time,  u(r,i,®)  is 
the  bandpass  filtered  synthetic  seismogram  with  central  frequency  ® ,  If;  (r'  -  r)  is  a  space 
window  centered  at  r ,  p  =  is  the  slowness  vector,  is  the  unit  vector  of  the  slowness 
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direction,  p^v  ‘  is  the  wave  slowness,  v  =  Vp  or  v  =  is  P-  or  ^-wave  velocity,  Q  is  a 
normalization  factor  determined  by  the  size  of  the  space  window  and  bandwidth  of  the 
frequency  filter,  etc.  The  space  window  is  for  a  small  seismic  array  with  size  related  to 
both  space  resolution  and  slowness  resolution.  A  larger  array  gives  better  slowness 
resolution  but  tends  to  smear  the  spatial  resolution,  while  a  small  array  gives  better 
spatial  resolution  but  less  accurate  slowness  calculation.  A  proper  trade-off  between 
space  and  slowness  resolution  is  required.  The  receiver  interval  should  be  small  enough 
to  avoid  spatial  aliasing.  The  average  energy  density  of  the  wavefield  as  a  function  of 
space,  time,  slowness  and  frequency  can  be  obtained  as 

D[r,ip,t,a))  =  pu^{r,-p,t,(o),  (2) 

and  the  energy  flux  related  to  slowness  vector  p  can  be  calculated  as 

J  =  vZt(r,p,?,£w)ep,  (3) 

where  p  is  the  density  and  m  (r,p,?,®)  is  the  rms  amplitude  of  the  stacked  velocity 
seismogram  u(r,p,^,®) .  Similarly,  we  can  slant  stack  seismograms  according  to  their 
horizontal  slowness  p;, 

u (r, p;,, ^ X (r*  -  r/,, z) u [r ' , z, ?  - p  Jr;  -  r;, ) ,  ®]  ,  (4) 

where  is  the  horizontal  position  and  r  =  (r;,,z) ,  is  a  2-dimensional  space  window, 
and  Cj  is  a  normalization  factor  similar  to  Q  in  equation  (1).  The  average  energy  density 
as  a  function  of  p,,  can  be  expressed  as 

^(r,P/,,^,®)  =  >OM^(r,p;,,^,®).  (5) 

The  pure  horizontal  energy  flux  related  to  horizontal  slowness  p;,  can  be  obtained 
from  equation  (3) 

Jjr,p;„^,®)  =  j  .  (6) 

In  the  waveguide,  the  energy  passing  through  a  surface  S  within  a  frequency  range  fl ,  a 
time  window  T  and  a  slowness  band  P  can  be  calculated  as 

E{s,P,T,a)  =  \ -tjwdtd^ds  ,  (7) 

where  e„  is  the  unit  normal  vector  of  the  surface  element  ds  ,  and  dp  =  dp^dpydp^ .  When 
choosing  a  vertical  intersection  as  the  surface,  equation  (7)  becomes 

E{S,P,T,PL)  =  \y.\y,\^\^J,dwdtdp,dz .  (8) 

Equations  (7)  and  (8)  provide  the  basis  for  extracting  energy  from  joint  domains. 
To  investigate  the  near-source  energy  partitioning  and  regional  phase  excitation,  the 
energy  generated  from  specific  mechanisms  is  decomposed  into  multiple  domains  and 
analyzed  based  on  its  dynamic  and  kinematic  properties.  The  condition  for  energy  to  be 
trapped  in  the  crustal  waveguide,  i.e.,  \^h\'^'^lvs-mantie  (where  is  the  upper  mantle  5 

velocity)  is  applied.  We  then  use  the  joint  window  {s,P,T,n)  to  sort  the  trapped  energy 
and  estimate  the  contribution  of  specific  mechanisms  to  regional  phases  such  as  Lg. 
Equations  (7)  or  (8)  can  be  partially  integrated,  which  allows  the  energy  to  be  projected 
onto  different  domains.  The  analysis  within  multiple  domains  provides  additional 
information  to  characterize  the  contributions  from  different  mechanisms.  This  method 
has  the  flexibility  that  we  can  either  intercept  the  entire  waveguide  energy  flux  or  just 
monitor  the  energy  from  specific  phases  or  mechanisms.  The  calculation  using  three- 
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dimensional  slowness  analysis  has  the  advantage  that  energy  is  fully  expanded  in  the 
entire  slowness  domain.  The  slowness  vector  p  for  both  P-  and  ^-waves  can  be  obtained 
independently,  giving  us  more  information  to  investigate  complicated  near-source 
processes.  The  calculation  based  on  horizontal  slowness  analysis  gives  the  energy  as  a 
function  of  horizontal  slowness  p;, .  It  provides  the  necessary  information  to  separate  the 
trapped  and  leaking  energy.  Although  it  does  not  directly  give  the  full  slowness  domain 
information,  by  combining  dynamic  and  kinematic  characteristics,  we  can  resolve  the 
near-source  phenomena  in  most  cases  without  ambiguity. 

Equations  (7)  and  (8)  are  expressed  with  energy,  since  energy  can  be  directly 
summed.  However,  regional  phase  observations  are  usually  taken  from  amplitudes. 

To  compare  the  numerical  prediction  with  observations,  we  calculate  the  normalized 
square  root  of  the  energy 

A^[EIE,f\  (9) 

where  E^  is  a  normalization  factor  which  can  be  seen  as  the  response  from  a  unit  source 
that  has  the  same  source  time  function  as  that  used  in  the  simulation.  Normalized  square 
root  energy  A  is  consistent  with  the  conventional  Zg-wave  measurement  based  on  the 
rms  amplitude  of  the  waveforms.  Another  advantage  of  using  A  is  that  the  variation  of 
amplitude  distribution  is  smoother  than  the  variation  of  energy.  For  these  reasons, 
throughout  this  report,  we  will  use  the  normalized  square  root  energy  in  all  figures 
although  we  sometimes  simply  call  it  “energy”. 
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Figure  1.  Different  slowness  domain  displays  for  2D  and  3D  models.  For  2D  model:  (a) 
energy  distribution  as  a  funetion  of  slowness  p,  and  (b)  as  a  funetion  of 
horizontal  slowness  Pf^  and  depth  z.  For  3D  model:  (e)  energy  distribution  as  a 
funetion  of  slowness  p  and  (d)  as  a  funetion  of  horizontal  slowness  p;,  and 
depth.  Energy  with  horizontal  slowness  larger  than  the  upper  mantle  S-wave 
slowness  will  be  trapped  in  the  waveguide  and  will  form  guided  S  phases  like  Lg. 

The  above  equations  are  for  3D  models.  In  a  2D  model,  the  y-component 
disappears  and  there  is  only  x  and  z  components  left.  Figure  1  shows  different  energy 
distributions  in  slowness  domains  for  2D  and  3D  models.  Depending  on  the  purpose,  the 
energy  can  either  be  expressed  as  a  function  in  slowness  domain,  or  as  a  function  in 
mixed  slowness-depth  domain.  The  upper  mantle  S-wave  slowness  Hvs-mantie  is  marked 
for  different  cases.  All  S'- wave  energy  with  horizontal  slowness  larger  than  the  upper 
mantle  S-slowness,  whether  directly  radiated  from  the  source  or  generated  as  secondary 
phases,  will  be  trapped  in  the  crustal  waveguide  and  will  contribute  to  the  guided  regional 
phases  unless  subsequent  scattering  causes  it  to  leak  out.  Figure  2  is  a  sketch  showing  the 
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configuration  of  the  local  slowness  analysis  method  for  investigating  the  near-source 
processes.  The  model  uses  an  explosion  source,  a  fine  scale  near-source  velocity  model 
and  a  short  distance  receiver  array  to  provide  synthetic  seismograms  for  the  ensuing 
slowness  analysis. 
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Figure  3.  Example  of  slowness  analysis  at  180  km  distance. 

Figure  3  gives  an  example  of  slowness  analysis.  Shown  on  the  top  is  the  synthetic 
seismogram  at  the  center  of  the  mini  array.  In  the  middle  are  energy  distributions  for  P- 
coda,  Lg  and  i?g-waves  in  2D  slowness  domain.  These  energy  distributions  clearly  show 
that  the  P-coda  is  composed  of  P-  and  reflected  /'S'-waves.  Both  of  these  have  relatively 
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small  horizontal  slowness.  The  Tg-wave  is  composed  of  multiply  reflected  5'- waves.  Its 
energy  falls  on  the  ^-wave  slowness  circle,  and  part  of  the  energy  stays  on  the  right  of  the 
upper  mantle  S-slowness  and  forms  trapped  phases.  The  i?g-wave  is  also  a  trapped  mode 
with  horizontal  slowness  larger  than  the  iS-slowness.  Shown  at  the  bottom  of  Figure  3  is 
the  energy  isolated  by  the  slowness  analysis.  The  horizontal  axis  is  time  or,  equivalently, 
the  inverse  of  the  group  velocity.  The  vertical  coordinate  is  the  horizontal  slowness  or, 
equivalently,  the  apparent  horizontal  phase  velocity.  The  solid  circles  are  energy 
measured  in  the  slowness  domain  with  their  sizes  being  proportional  to  the  amount  of 
energy.  The  horizontal  dashed  line  marks  the  upper  mantle  S-slowness  that  divides  the 
trapped  and  leaky  energy. 


(a)  0.3-1. 5  Hz 
group  velocity  km/s 

5,0  4.0  3.5  3.0 


(b)  2.0-5.0  Hz 
group  velocity  km/s 


time  sec  time  sec 

Figure  4.  Slowness  analysis  calculated  for  the  EK  model  and  EK  model  with  random 
patches  at  different  depths.  Column  (a)  is  for  frequency  band  0.3- 1.5  Flz,  (b)  is 
for  2. 0-5.0  Flz.  The  top  row  is  for  background  velocity  model  and  the  lower  rows 
are  for  models  with  random  patches.  The  depths  of  random  patches  are  labeled  in 
the  panels.  Details  see  the  text. 


3.3.  The  Energy  Partitioning  Formalism 

For  convenience,  we  symbolically  write  the  near-source  energy  partitioning  process 
for  an  explosion  source  as 

E’^{p,f)  =  S{f)R'^{p,f),  (10) 

where  {p,f)  is  the  near-source  energy  partitioned  to  the  type  K  wave  {K  can  be  P,  S, 
Lg  ,  Rg  or  other  wave  types),  p  is  the  slowness,  /  is  the  frequency,  S{f)  is  the  spectrum 
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of  an  isotropic  explosion  source.  [p,f)  is  the  energy  response  funetion  of  the  near 
souree  strueture  for  exeiting  type  K  wave,  and  ean  be  expressed  as 

R^{pJ)  =  R’^F{pJ)  +  I.4{pjy^\pJ)-T.RF{pJ)T'^"\pJ)-  (11) 

j  .1 

On  the  right  hand  side  of  this  equation,  Rf  [p,f)  is  the  response  of  a  flat,  homogeneous 
layered  earth  model,  partitioning  the  souree  energy  into  different  phases.  The  transfer 
funetion  provides  the  J-to-K  eoupling  which  modifies  the  original  partitioning  by 
moving  energy  from  one  phase  to  another.  The  second  term  on  the  right  hand  side 
denotes  energy  being  imparted  into  the  K  wave  through  eoupling  and  the  third  term 
denotes  energy  lost  from  the  K  wave  to  other  phases.  The  eombined  effeet  gives  the  total 
partitioning  of  the  energy  radiated  from  an  isotropie  souree  into  the  K  wave  energy 
distributed  in  slowness  and  frequeney  domains.  This  energy  will  develop  into  different 
regional  phases  which  propagate  to  remote  distanees.  Having  a  eomplete  deseription  of 
the  slowness  distribution  allows  us  to  aeeurately  prediet  energy  imparted  to  the  distant 
regional  phases  based  on  the  near-souree  energy  budget. 

Table  1.  Eastern  Kazakh  Model 


top  of  layer  (km) 

Vp  (km/s) 

Vs  (km/s) 

p  (gram/em  ) 

0.0 

5.05 

2.91 

2.70 

1.0 

5.66 

3.21 

2.70 

5.0 

5.55 

3.17 

2.70 

10.0 

6.01 

3.41 

2.70 

15.0 

6.29 

3.55 

2.90 

25.0 

6.65 

3.73 

2.90 

30.0 

6.89 

3.85 

2.90 

45.0 

7.59 

4.20 

3.00 

47.5 

8.33 

4.57 

3.30 

Investigating  these  response  and  transfer  funetions  provides  a  way  to  estimate  the 
underlying  proeess  of  energy  partitioning.  Due  to  the  diverse  meehanisms  involved,  the 
aetual  near-souree  energy  partitioning  ean  be  highly  eomplex.  Faetors  sueh  as  the  souree 
depth,  loeal  layered  strueture,  attenuation,  random  volumetrie  veloeity  perturbations  and 
free  surfaee  fiuetuations  all  affect  the  partitioning.  These  effects  often  contribute  to  the 
partitioning  in  a  coupled  way  and  the  entire  proeess  is  not  neeessarily  linear  or  simply 
separable.  In  this  study,  we  foeus  on  the  eontribution  from  surfaee  seattering.  We  use 
aeeurate  numerieal  modeling  to  simulate  the  eomplex  partitioning  proeess  and  the 
slowness  analysis  to  ealeulate  the  diserete  response  funetions.  The  above  symbolic 
equations  provide  us  with  a  basie  formalism  for  understanding  the  proeess.  The 
embedded  array  method  eannot  eharaeterize  all  partitioning  eoeffieients,  but  is 
partieularly  well  suited  to  eharaeterizing  relative  ehanges  in  the  trapped  energy 
distribution  for  waves  that  will  travel  to  large  distanees  in  the  crustal  waveguide. 
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waveguide  energy  surface  energy 

flux  at  1 80  km  at  450  km 


0.0  2.5  5.0  7.5  10.0  12.5  0.0  2.5  5.0  7.5  10.0  125 


0.0  2.5  5.0  7.5  10.0  12.5  0.0  2.5  5.0  7.5  10.0  12.5 

depth  of  randomness  km  depth  of  randomness  km 


Figure  5.  Comparison  between  waveguide  energy  flux  at  180  km  (left  column)  and  the 
wave  energy  on  the  surface  at  450  km  (right  column)  for  Pg,  Lg  and  Rg 
windows.  The  frequency  range  is  0.3- 1.2  Hz.  Shown  in  each  panel  is  relative 
energy  versus  the  depth  of  random  patches.  Dashed  lines  indicate  the  energy 
level  for  background  model  and  short  bars  indicate  the  energy  changes  due  to  the 
near-source  scattering.  If  the  FDSA  method  works,  the  two  columns  should  be 
consistent. 

3.4.  Testing  the  Validity  of  the  Method 

We  will  first  check  the  validity  of  the  FDSA  method  by  comparing  the  energy 
partitioning  predicted  at  short  distances  using  slowness  analysis  with  the  energy 
measured  from  long  distance  surface  receivers  using  a  conventional  method.  For  all 
numerical  examples  calculated  in  this  paper,  unless  otherwise  indicated,  we  use  the 
horizontally  layered  Eastern  Kazakh  (EK)  model  (Priestley  et  ah,  1988)  as  the 
background  and  modify  it  by  adding  random  velocity  fluctuations  at  different  locations. 
The  EK  model  (Table  1)  has  a  high  velocity  top  layer  withv^  =  5.05  km/s  and  an  upper 
mantle  S'-wave  velocity  vs^^amie  =  4.57  km/s.  Eor  this  layered  model,  E-waves  radiated 
from  an  explosive  source  cannot  be  effectively  converted  to  trapped  ^-waves.  The 
random  velocity  perturbations  added  have  an  exponential  power  spectrum  with  horizontal 
and  vertical  correlation  lengths  equal  to  0.5  km.  We  call  the  part  of  the  model  with 
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velocity  perturbations  the  “random  patch”.  Within  the  patch,  both  P-  and  iS-wave 
velocities  have  the  same  relative  rms  perturbation  and  the  relative  density  perturbation  is 
50%  of  the  relative  velocity  perturbation.  To  eliminate  the  possible  effect  of  sharp  edges 
of  these  random  patches,  a  space  window  with  smoothed  edges  is  applied  to  the  patch. 
The  sizes,  locations  and  the  rms  perturbations  of  these  patches  will  vary  depending  on  the 
purpose  of  the  investigation. 

We  first  test  models  with  scattering  at  different  depths  by  varying  the  location  of 
random  patches  within  the  EK  model.  These  patches  are  2.5  km  in  vertical  extent,  extend 
horizontally  from  5  to  25  km  and  are  located  at  different  depths.  Within  the  patches,  the 
P-  and  ^-wave  velocities  have  10%  rms  fluctuations.  An  explosion  source  is  located  at 
depth  0.5  km.  Figure  4  shows  the  slowness  analysis  results  at  a  distance  of  180  km  in  the 
crustal  waveguide.  Columns  (a)  and  (b)  are  for  frequencies  0.3-1. 5  Hz  and  2. 0-5.0  Hz, 
respectively.  Each  panel  is  similar  to  the  example  shown  in  Figure  3.  The  top  row  is  for 
the  background  velocity  model  and  the  other  rows  are  calculated  for  models  with  near¬ 
source  random  velocity  perturbations.  The  depths  of  the  random  patches  are  labeled  in 
the  figure.  We  first  focus  on  the  low-frequency  results  in  column  (a).  For  the  background 
velocity  model,  there  is  a  strong  Rg  phase  but  very  little  energy  within  the  Lg  group 
velocity  window  of  3. 0-3. 5  km/s,  which  is  typical  for  a  crustal  model  with  a  high  velocity 
top  layer.  For  models  with  random  velocity  patches,  compared  with  the  background 
model,  considerable  energy  is  transferred  to  the  Fg-wave  through  scattering,  while  the 
i?g-wave  is  weakened.  The  tendency  is  that  the  shallower  the  random  patch,  the  more 
energy  is  scattered  into  Lg.  For  the  EK  model,  which  has  a  thick  high  speed  crust,  the 
distance  for  Fn  to  cross  Pg  at  the  free  surface  is  about  200  km.  At  180  km  distance,  there 
is  no  prominent  Pn  energy  shown  in  this  analysis. 

For  high  frequency  results  in  column  (b),  from  the  top  panel  we  see  strong  Pg 
energy.  There  is  also  energy  within  the  Lg  group  velocity  window.  However,  this  energy 
has  a  similar  horizontal  slowness  to  the  Fg-wave,  implying  that  it  is  generated  from  Pg 
through  P-io-P  and  P-io-S  reflections  on  the  free  surface  and  interfaces  such  as  the  Moho 
discontinuity.  Although  the  energy  exists  at  short  distances,  the  steep  incident  angle 
causes  energy  to  gradually  leak  to  the  upper  mantle  through  multiple  reflections  and  it 
cannot  form  trapped  regional  phases.  In  the  other  panels,  after  adding  random  velocity 
perturbations  in  the  near-source  region,  part  of  the  P  energy  transfers  to  the  Lg  wave,  i.e., 
energy  falls  into  the  proper  group  velocity  and  slowness  windows.  In  general,  the 
scattering  affects  Rg-io-Lg  energy  at  low  frequencies  and  affects  P-io-Lg  coupling  at  high 
frequencies.  Summing  up  energy  through  the  waveguide  cross  section  and  within  the 
proper  time  (group  velocity)  windows  allows  the  energy  for  the  related  wave  types  to  be 
obtained.  Specifically,  summing  up  the  energy  located  within  the  Lg  group  velocity 
window  and  above  the  upper  mantle  S-slowness  allows  Lg  energy  at  long  distances  to  be 
predicted.  It  is  necessary  to  confirm  this  calculation. 

Figure  5  compares  the  waveguide  energy  obtained  using  different  methods.  The 
left  panel  shows  the  square  root  energy  passing  through  the  waveguide  cross  section  at 
180  km  calculated  from  the  slowness  analysis  shown  in  Figure  4  and  the  right  panel 
shows  the  rms  amplitudes  at  450  km  calculated  using  a  conventional  processing 
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technique  (i.e.,  computing  surface  synthetic  seismograms,  applying  frequency  filters  and 
group  velocity  windows  to  isolate  different  phases,  and  calculating  rms  amplitudes  for 
these  phases.).  The  same  velocity  model  used  in  Figure  4  is  adopted  here  and  a  full-scale 
finite-difference  simulation  is  computed  to  provide  data  up  to  500  km  distance.  For 
simplicity,  we  label  the  predicted  trapped  energy  as  “Lg”  energy.  The  relative  energy 
changes  of  Pg-,  Lg-  and  Rg-waves,  as  functions  of  the  depth  of  the  random  velocity 
patches  are  shown  in  Figure  5.  The  calculations  for  different  phases  demonstrate  that  the 
relative  change  of  energy  obtained  in  the  waveguide  slowness  analysis  corresponds 
closely  to  that  obtained  on  the  free  surface,  even  for  dramatic  changes  such  as  the  strong 
scattering  of  Rg.  The  results  confirm  that  the  slowness  analysis  within  the  waveguide 
correctly  predicts  the  surface  regional  observations  at  greater  distance. 


Figure  6.  Examples  showing  the  slowness  analyses  at  distances  (a)  50  km  and  (b)lOO  km. 

The  and  are  upper-mantle  S-slowness,  crustal  P-slowness  and  S- 

slowness,  respectively.  The  energy  that  can  be  trapped  in  the  crustal  waveguide  is 
indicated  in  the  figure. 
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Figure  7.  Comparison  between  trapped  waveguide  energy  measured  at  50  km  (horizontal 
coordinate)  and  100  km  (vertical  coordinate)  for  frequency  bands  0.5-1. 5  Flz, 
1. 0-3.0  Flz  and  2. 0-4.0  Flz.  Different  dots  are  results  from  different  velocity 
models  and  source  depths.  The  results  show  a  general  linear  relationship  for  all 
frequencies. 
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To  examine  the  energy  flux  measurements  obtained  at  even  shorter  distanees,  we 
eompare  the  slowness  analysis  measurements  at  50  km  with  those  at  100  km.  For  these 
two  distanees,  Figure  6  gives  examples  of  horizontal  slowness  analyses,  whieh  show 
quite  different  features  in  the  slowness-depth  domain  due  to  evolution  of  the  wavefield 
with  range.  Figure  7  eompares  the  eorresponding  waveguide  energy  measured  at  these 
two  distanees  for  different  frequeney  bands.  The  vertieal  and  horizontal  eoordinates  are 
for  measurements  at  the  two  distanees  and  dots  represent  results  for  different  souree 
depths  and  veloeity  models.  Although  a  wide  range  of  near-souree  struetures  and  souree 
depths  are  used  to  generate  these  measurements,  the  results  show  a  general  linear 
relationship  for  all  frequeney  bands.  This  further  verifies  that  we  ean  use  a  small  model  to 
investigate  the  near-souree  energy  partitioning  robustly. 

4.  Regional  Phase  Excitation:  Contributions  from  Volumetric  Scattering 

To  investigate  seattering  eontributions  from  volumetric  heterogeneities,  we  use  a 
2D  finite-difference  method  to  simulate  the  wave  propagation.  The  model  uses  an 
explosion  source,  a  fine  scale  near-souree  velocity  model  and  a  short  distance  receiver 
array  to  provide  synthetic  seismograms  for  the  ensuing  slowness  analysis  (see  Figure  2). 

4.1.  P-pS-to-Lg  and  P-to-Lg  Conversion 

We  first  investigate  the  P-pS-io-Lg  conversion  caused  by  near-souree  lateral 
velocity  variations  and  assess  its  effect  on  the  explosion  ^-wave  energy  budget.  In  a 
horizontally  layered  model  with  overburden  P-wave  velocity  larger  than  the  upper  mantle 
iS-wave  velocity,  the  free  surface  reflected  j^^-wave  has  a  steep  incidence  angle  and 
cannot  be  trapped  in  the  crustal  waveguide  to  form  Lg.  In  this  case,  the  energy  transfer 
through  P-pS-io-Lg  coupling  is  almost  zero.  Although  it  is  generally  agreed  that  the 
existence  of  near-souree  lateral  velocity  variation  can  increase  the  P-io-Lg  energy 
exchange,  the  detailed  mechanism  underlying  this  process  is  still  not  fully  understood. 

Figure  8  compares  the  simulated  P-pS-io-Lg  coupling  in  models  with  and  without 
near  surface  lateral  velocity  variations.  Figure  8a  is  for  the  EK-model.  A  shallow 
explosion  source  located  at  depth  0.5  km  generates  P,  pS  and  Rg  waves.  Two- 
dimensional  slowness  analysis  is  conducted  for  selected  phases  in  the  wavefield  and  the 
results  are  shown  together  with  the  wavefield  snapshot.  The  synthetic  seismograms  were 
bandpass  filtered  between  2.0  and  6.0  Hz  before  the  slowness  analysis.  As  can  be  seen 
from  the  result,  the  P-wave  leads  the  wavefield  and  has  a  distinct  slowness. 
Reverberations  within  the  uppermost  crust  causes  multiply  parallel  pS  wavefronts  with 
their  horizontal  slowness  approximately  equal  to  the  overburden  R-slowness.  The  pS 
energy  stays  to  the  left  of  the  upper-mantle  ^-slowness  and  there  is  no  energy  transferred 
from  P  to  Lg.  In  Figure  8b,  a  shallow  random  velocity  patch  is  added  to  the  EK-model  to 
test  the  effect  of  near  surface  scattering.  The  random  patch  has  a  5%  rms  velocity 
fluctuation  and  is  located  at  distances  5  to  15  km  and  depths  0  to  2.5  km  (shown  in  the 
snapshot  as  a  shaded  area).  The  slowness  analyses  are  conducted  fov  P,pS  and pS-coda. 
Although  P-  and  early  pS-waves  are  barely  affected,  the  pS-coda  clearly  contains  some 
scattered  energy  with  horizontal  slowness  to  the  right  of  the  upper-mantle  ^-slowness. 
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To  investigate  further  the  seattering  from  a  shallow  random  pateh,  horizontal 
slowness  analysis  is  eondueted  at  a  distanee  of  20  km  and  depth  0-12.5  km.  Figure  9 
shows  the  energy  distribution  in  the  slowness-depth  domain  with  arrival  times  and  major 
phases  labeled  in  the  frames.  The  two  prominent  down-going  phases  are  the  P-wave  and 
the  free-surfaee  refleeted pS-wsive.  The  Rg  energy  enters  the  array  at  6.0  s  with  its  depth 
close  to  the  surface  and  slowness  beyond  the  ^-slowness.  Shown  in  row  (a)  is  the  result 
using  the  EK-model.  Due  to  the  nearly  horizontal  propagation  of  the  P-wave  at  the  free 
surface,  the  pS-wa\e  energy  has  a  horizontal  slowness  which  is  similar  to  the  overburden 
P-slowness  and  the  energy  falls  to  the  left  of  the  upper  mantle  S'-slowness.  In  rows  (b) 
and  (c),  shallow  random  velocity  patches  with  rms  velocity  fluctuations  3%  and  5%  are 
added  to  the  EK-model  at  distances  5  to  15  km  and  depths  0  to  2.5  km  (the  same  position 
as  in  Eigure  8b).  As  can  be  seen  in  the  figure,  scattering  causes  part  of  the  pS  energy  to 
cross  the  upper  mantle  ^-slowness  and  build  up  in  the  dashed  rectangles.  At  shallow 
depths,  the  slowness  of  the  scattered  energy  approaches  the  ^-slowness.  With  an  increase 
in  depth,  the  slowness  of  this  energy  gradually  merges  with  the  E-slowness.  The  slowness 
behavior  is  consistent  with  scattering  of  waves  at  shallow  depths,  increasing  as  shallow 
heterogeneity  increases. 


(a)  EK-model 

pS  coda  pS  P 


(b)  EK+random 

pS  coda  pS  P 


Figure  8.  P-pS-Lg  conversion  due  to  shallow  scattering,  (a)  is  for  the  EK-model  and  (b)  is 
for  the  EK-model  with  a  shallow  random  patch.  The  random  patch  has  a  5%  rms 
velocity  fluctuation  and  is  shown  as  a  shaded  area  between  horizontal  distances  5 
to  15  km  and  depths  0  to  2.5  km.  The  snapshots  and  results  of  slowness  analyses 
of  P,  pS  and  pS-coda  are  shown  in  the  figure.  Details  are  given  in  the  text. 


Eigure  10  investigates  scattering  taking  place  at  deeper  depths.  The  configuration 
of  the  source  and  model  is  similar  to  that  used  in  Eigure  8b,  except  the  random  patch  with 
3%  rms  velocity  fluctuation  is  added  to  the  EK-model  between  distances  5  and  15  km 
and  depths  2.5  and  10  km  (shown  in  the  snapshot  as  a  shaded  area).  The  2D  slowness 
analysis  is  conducted  for  selected  phases  in  the  wavefield  and  the  results  are  presented  in 
the  figure.  After  passing  through  the  random  region,  there  is  E-coda  composed  of 
scattered  E-  and  ^-waves  generated  from  the  direct  E-wave.  Although  the  early  part  of 
the  pS-wavQ  does  not  contribute  to  the  trapped  energy,  its  later  part  contains  energy 
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located  to  the  right  of  the  upper-mantle  S-slowness  whieh  therefore  will  contribute  to  the 
trapped  regional  phases.  Figure  1 1  gives  the  energy  distribution  in  the  slowness-depth 
domain  for  different  models  where  row  (a)  is  for  the  EK-model  and  rows  (b)  and  (e)  are 
for  the  EK-model  with  3%  and  5%  rms  fluetuations  in  a  random  patch  like  that  used  in 
Eigure  10.  The  slowness  analysis  is  eondueted  at  a  distance  of  20  km  and  for  depths 
between  0  and  12.5  km.  As  expeeted,  with  the  EK-model  no  energy  is  seen  beyond  the 
upper  mantle  ^-slowness,  but  after  the  lateral  veloeity  variations  are  introdueed,  energy 
starts  to  build  up  to  the  right  of  the  upper-mantle  S'-slowness.  Two  types  of  seattered 
energy  ean  be  found  in  the  slowness-depth  domain:  weak  but  widely  distributed  S  energy 
(indieated  by  the  dashed  ellipses)  and  seattered  energy  linked  to  the  pS-wave  (indieated 
by  the  dashed  reetangles).  Both  types  of  energy  satisfy  the  eriterion  >  Ps-mamie  will 
eontribute  to  the  Eg-wave.  The  widely  spread  scattered  S-wavc  is  generated  by  the  E-to- 
Lg  eoupling  through  volumetric  scattering.  The  seattering  proeess  redistributes  the  angle 
speetrum  of  the  original  ineident  waves.  Sinee  there  is  no  lateral  heterogeneity  at  the  top 
of  the  crust,  scattered  pS-wavQS  are  generated  by  the  interaction  between  distorted 
ineident  E-waves  and  a  smooth  free  surface. 

Comparing  Eigure  1 1  to  Eigure  9,  the  seattered pS-wavQ  generated  from  a  deeper 
random  pateh  appears  later  than  that  from  a  shallow  random  pateh.  Both  volumetrie 
scattering  and  seattering  near  the  free-surface  affect  the  general  P-to-Lg  conversion. 
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Figure  9.  Florizontal  slowness  analyses  for  investigating  the  P-pS-Lg  coupling  with  (a) 
the  EK-model,  (b)  the  EK-model  plus  a  3%  shallow  random  patch  and  (c)  the 
EK-model  plus  a  5%  random  patch.  The  slowness  analyses  are  conducted  at 
distance  20  km  and  for  depth  range  0  to  12.5  km.  The  configuration  of  the  source 
and  model  is  the  same  as  that  used  in  Figure  8.  Major  phases  are  labeled  in  the 
figure  and  energy  circled  by  dashed  rectangles  is  scattered  pS-wave. 


15 


depth  km 


pS  coda  pS  l3tei-  p-coda  early  P-coda 


distance  km 

Figure  10.  The  P-pS-Lg  conversion  due  to  a  deeper  random  patch.  The  random  patch  has 
a  3%  rms  velocity  fluctuation  and  is  shown  as  a  shaded  region  between  distances 
5  and  15  km  and  depths  2.5  and  10  km.  The  slowness  analysis  for  P,  P-coda,  pS 
and  pS-coda  are  also  shown  in  the  figure. 
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Figure  11.  Energy  distribution  in  depth  and  horizontal  slowness  domain  for  (a)  the  EK- 
model,  (b)  EK-model  plus  a  3%  random  patch  and  (c)  EK-model  plus  a  5% 
random  patch.  The  configuration  of  the  source  and  model  is  same  as  that  used  in 
Figure  10.  Energy  circled  by  dashed  rectangles  is  P-pS-Eg  scattering  and  energy 
circled  by  dashed  ellipses  is  P-Eg  scattering. 


16 


0  10  20  30  40  SO  60 


distance  km 


(b)  source  depth  3.0  km 


Figure  12.  Wavefield  snapshots  for  explosion  sources  at  depths  (a)  0.5  km  and  (b)  3.0  km.  Note 
that  a  shallower  explosion  is  a  more 

4.2.  Contributions  from  the  S*-Wa\e 

For  shallow  explosion  sources,  the  5'*-wave  may  become  a  significant  contributor 
to  Lg  (Gutowski,  et  al.,  1984;  Lilwall,  1988;  Xie  and  Lay,  1994;  Vogljord,  1997).  The 
amplitude  of  S'*  can  be  large  if  the  source  depth  is  within  a  fraction  of  a  wavelength  from 
an  interface.  This  makes  its  excitation  highly  dependent  on  the  source  depth  and 
frequency.  Figure  12  shows  snapshots  for  explosion  sources  at  0.5  km  and  3.0  km, 
respectively.  The  result  clearly  shows  that  a  shallow  source  generates  larger  S*-  and  Rg- 
waves.  We  investigate  the  contribution  of  the  S*-wave  within  the  EK  model.  Figure  13 
shows  horizontal  slowness  analyses  at  a  distance  35  km  and  for  depths  between  0  and  30 
km.  The  time  window  is  chosen  between  1 1  to  13  s  after  the  direct  E-wave  passes  the 
receiver  array.  The  synthetic  seismograms  are  bandpass  filtered  between  1.0  -  5.0  Hz. 

The  four  rows  from  top  to  the  bottom  correspond  to  source  depths  0.25  km,  0.5  km,  1.0 
km,  and  2.0  km,  respectively.  The  major  arrival  is  the  down-going  free  surface  reflected 
pS-wawe,  which  has  a  horizontal  slowness  similar  to  the  overburden  E-slowness.  As 
expected  for  a  horizontally  layered  model,  the  pS  energy  stays  to  the  left  of  the  upper 
mantle  N-slowness  and  has  no  contribution  to  the  trapped  regional  phases.  For  shallow 
sources,  the  Eg-wave  enters  the  array  at  about  12  s  and  its  energy  concentrates  between  0 
to  3  km,  as  can  be  seen  on  the  upper  right  corners  in  the  slowness-depth  domain.  For 
source  depth  of  2.0  km,  the  Eg-wave  is  very  weak.  The  5'*-wave  enters  the  array  from  a 
shallow  depth  and  gradually  merges  with  the  pS-wavQ  (also  refer  to  Figure  12).  The  S*- 
wave  is  strong  for  shallow  sources  and  its  amplitude  decreases  with  increasing  source 
depth.  Very  little  S*  energy  can  be  observed  for  source  depths  below  2  km.  In  the  joint 
domains,  the  5'*-energy  can  be  isolated  and  quantified  even  within  a  complicated 
wavefield,  which  is  very  difficult  using  remote  surface  synthetics.  The  dashed  rectangles 
are  the  time-slowness-depth  window  used  to  locate  the  S*  energy.  The  time  window  is 
chosen  after  the  arrival  of  direct  E-wave,  and  a  variable  depth  range  is  chosen  to  avoid 


17 


contamination  from  the  Rg-wave.  The  slowness  range  is  chosen  between  0.23  and  0.34 
sikm  .  The  energy  from  successive  windows  ean  be  summed  together  to  give  the 
eontribution  of  S'*  to  the  trapped  regional  phases. 


efficient  source  for  generating  S*  and  Rg-waves. 


source 

depth 

0.25  km 

0.50  km 

1.00  km 


Figure  13.  Slowness  analysis  for  investigating  S*-to-Lg  conversion.  Different  rows  are 
for  different  source  depths.  Dashed  rectangles  indicate  the  time-space-slowness 
windows  used  to  pick  the  S*  energy. 


4.3.  The  Frequency  Dependent  Lg  Excitation  Function 

The  frequency  dependence  of  Tg-wave  exeitation  is  rooted  in  the  underlying 
physieal  processes  and  is  usually  eontrolled  by  different  characteristic  scales.  For 
example,  the  excitation  of  Lg  by  5**,  Rg-io-S  scattering  and  spall  are  all  highly  souree 
depth  dependent.  The  excitation  spectra  from  individual  or  joint  mechanisms  contributing 
to  regional  phases  depiet  the  frequency  dependence  of  these  processes.  Frequeney 
dependent  PIS  ratios  will  depend  on  the  excitation  functions  of  multiple  phases.  We  use 
FDSA  to  quantify  Tg-wave  excitation  spectra  from  S*-  and pS-v^avQS.  Figure  14  gives 
the  S'*-to-Tg  exeitation  speetra  as  functions  of  source  depth  and  frequency.  The  slowness 
analysis  and  multi-domain  window  used  to  pick  the  trapped  energy  are  similar  to  that 
shown  in  Figure  13.  A  series  of  bandpass  fdters  is  used  to  give  responses  at  different 
frequeneies.  The  vertical  coordinate  is  the  normalized  relative  energy  .  Since  the 

source  time  function  has  been  taken  away,  the  exeitation  function  is  the  impulse  response 
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of  the  model  to  the  souree.  The  results  clearly  show  that  the  S*-to-Lg  excitation  is 
generally  enhanced  for  lower  frequency  and  shallow  source  depth.  The  major 
contribution  comes  from  sources  located  above  1  km.  For  sources  at  depths  below  1  km, 
only  low  frequency  energy  below  1  Hz  has  significant  contribution  to  Tg-wave 
excitation.  However,  the  responses  are  also  model  dependent.  For  a  model  with  a 
homogeneous  crust  (Figure  14a),  the  distribution  has  simple  monotonic  tendencies  in 
both  source  depth  and  frequency.  For  the  EK  model  (Figure  14b),  the  excitation  spectrum 
has  a  maximum  at  depth  1  km  and  a  more  complicated  frequency  dependence.  This  may 
reflect  the  fact  that  the  EK  model  has  an  interface  at  1  km  depth.  The  S'* -waves 
generated  or  reflected  from  multiple  interfaces  may  interfere  with  each  other  and  give  a 
complicated  frequency  spectrum. 


(a)  Constant  crustal  model  (b)  EK-model 


Figure  14.  Normalized  Eg  excitation  spectra  for  sources  in  different  velocity  models  and 
at  different  depths  with  (a)  a  model  with  a  homogeneous  crust  and  (b)  the  EK- 
model. 

To  investigate  the  combined  effect  for  S*-wave  and  near-source  scattering,  we 
add  shallow  random  velocity  patches  to  the  EK  model.  The  random  patch  extends 
between  distances  of  5  to  25  km  and  depths  of  0  to  2.5  km.  Shown  in  Eigures  15a  and  b 
are  excitation  spectra  for  random  patches  with  rms  velocity  fluctuations  of  3%  and  5%, 
respectively.  The  most  prominent  feature  is  the  build  up  of  high  frequency  energy.  The 
scattered  energy  increases  with  rms  velocity  fluctuations.  Eigures  15c  and  d  isolate  the 
scattered  energy  by  subtracting  the  excitation  spectrum  of  the  EK-model  from  the  spectra 
for  models  with  random  velocity  patches.  Two  types  of  energy  can  be  identified  within 
the  frequency-depth  domain.  The  high-frequency  energy  results  from  P-pS-to-Lg  and  P- 
to-Lg  scattering.  This  energy  is  especially  important  for  deeper  sources  to  generate  Lg- 
waves,  since  a  deeper  source  generates  little  trapped  energy  in  a  horizontally  layered 
model.  The  low-frequency  energy  concentrated  at  shallow  source  depths  comes  from  Rg- 
to-Lg  scattering. 
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(a)  EK+3%  shallow  random  (b)  EK+5%  shallow  random 


(c)  EK+3%  scattered  energy  (d)  EK+5%  scattered  energy 


Figure  15.  Normalized  Lg  excitation  spectra  for  sources  in  the  EK-model  with  shallow 
random  patches,  for  (a)  the  EK-model  with  a  3%  shallow  random  patch,  (b)  the 
EK-model  with  a  5%  shallow  random  patch,  (c)  and  (d)  the  isolated  scattered 
energy  in  (a)  and  (b)  respectively  due  to  the  random  patches  found  by  removing 
the  energy  for  the  layered  models.  Note  different  vertical  scales  are  used  for 
scattered  energy. 
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(c)  EK+3%  scattered  energy  (d)  EK+0.05  scattered  energy 


Figure  16.  Normalized  Lg  excitation  spectra  for  sources  in  the  EK-model  with  deep 
random  patches,  for  (a)  the  EK-model  with  a  3%  deep  random  patch,  (b)  EK- 
model  with  a  5%  deep  random  patch,  (c)  and  (d)  the  isolated  scattered  energy  in 
(a)  and  (b)  due  to  the  random  patches. 


Figure  16  gives  the  excitation  spectra  for  the  EK-model  with  deeper  random 
patches.  The  random  patch  is  located  between  distances  5  to  25  km  and  depths  7.5  to  10.0 
km.  Figures  16a  and  b  give  excitation  spectra  for  random  patches  with  rms  velocity 
fluctuations  of  3%  and  5%,  respectively.  Figures  16c  and  d  give  the  isolated  scattered 
energy.  The  scattered  energy  from  the  deeper  random  patches  has  little  low-frequency 
content,  which  supports  the  interpretation  that  the  low-frequency  energy  comes  from  the 
Rg-io-Lg  scattering.  The  frequency  dependent  excitation  spectra  establish  the  relationship 
between  the  observations  and  the  characteristics  of  sources  and  near-source  structures. 
They  provide  the  basis  for  evaluating  the  dominant  mechanisms  for  Tg-wave  excitation. 
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5.  Regional  Phase  Excitation:  Contributions  from  Surface  Scattering 

Figure  17  shows  the  basie  eonfiguration  of  the  ealculation.  Random  topographic 
fluctuations  described  by  an  exponential  power  spectrum,  a  correlation  distance  and 
prescribed  rms  fluctuation  are  used  for  the  free  surface.  We  use  a  2D  P-SV  wave 
boundary-element  method  (Ge  et  al.,  2005)  to  generate  synthetics  seismograms  for  an 
isotropic  explosion  source.  The  boundary  element  method  has  been  proven  to  be  an 
accurate  wave  propagation  approach  when  the  earth  model  includes  irregular  surfaces. 
The  synthetic  seismograms  are  computed  for  an  embedded  vertical  array  of  41x90 
receivers  located  between  distances  30  and  50  km  and  depths  0  and  45  km  (see  Figure  1). 
The  synthetic  seismograms  from  the  array  are  processed  using  the  local  slowness  analysis 
method. 
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Figure  17.  Configuration  of  the  source,  model  and  receiver  array. 


5. 1.  Phenomena  Resulting  from  Free  Surface  Scattering 

Space  domain  representation  of  surface  scattering 

Illustrated  in  Figure  18  are  boundary  element  generated  snapshots  for  wavefields  in  a 
model  with  free  surface  scattering.  The  parameters  of  the  three-layer  velocity  model  are 
listed  in  Table  2.  A  random  fluctuation  with  a  correlation  length  of  0.5  km  and  a  rms 
fluctuation  of  0.15  km  is  used  for  the  free  surface  and  extends  between  distances  30  and 
50  km.  The  source  is  located  at  distance  20  km  and  a  depth  of  0.5  km.  Figures  18a  and 
18b  are  horizontal  and  vertical  displacements.  In  addition  to  familiar  major  phases  (e.g., 
P,  pS,  Rg)  expected  for  a  shallow  explosion  in  a  flat  earth  model,  scattered  body  and 
surface  waves  from  the  rough  free  surface  are  present  in  the  wavefield.  The  surface-to- 
body  and  body-to-body  wave  scattering  is  distributed  through  the  entire  medium 
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following  the  direct  waves,  and  the  body-to-surface  and  surface-to-surface  wave 
scattering  is  concentrated  at  very  shallow  depths  following  direct  waves  as  they  graze  the 
surface. 

Table  2.  Three-layer  velocity  model 


Bottom  of  layer  Fp 

Ts 

P 

(km)  (km/s) 

(km/s) 

(g/cm^) 

10  5.6 

3.2 

2.7 

45  6.5 

3.8 

2.9 

Infinity  8.0 

4.5 

3.3 

Table  3.  Two-layer  velocity  model 

Bottom  of  layer  Fp 

Ts 

P 

(km)  (km/s) 

(km/s) 

(g/cm^) 

45  6.5 

3.6 

2.9 

infinity  8.0 

4.5 

3.3 

Figure  18.  Wavefield  snapshot  at  t  =  10.0  s  for  a  model  with  random  free  surface 
fluctuation,  where  (a)  and  (b)  are  horizontal  and  vertical  components  of  the 
wavefield,  and  (c)  and  (d)  are  horizontal  and  vertical  components  of  the  scattered 
wavefield  obtained  by  subtracting  the  flat  model  wavefield  from  the  random 
surface  wavefield.  The  source  is  at  20  km  horizontal  position  and  0.5  km  depth. 

To  isolate  the  scattered  phases,  we  subtract  the  wavefield  generated  for  a  flat  surface 
from  the  wavefield  for  the  model  having  a  rough  surface,  yielding  the  results  presented  in 
Figures  18c  and  18d.  Most  of  the  scattered  body  waves  are  shear  waves.  Due  to  the 
coupling  between  different  wavenumbers,  the  scattered  body  waves  have  a  very  broad 
range  of  propagation  directions.  The  horizontal  component  includes  mostly  shear  waves 
propagating  with  steep  dip  angles  that  will  tend  to  contribute  to  teleseismic  ^-waves.  The 
shear  waves  on  the  vertical  component  mostly  have  shallow  angles  and  will  contribute  to 
crustal  guided  waves  such  as  regional  Lg.  Although  these  space  domain  images  are 
instructive  for  understanding  some  aspects  of  the  surface  scattering,  to  fully  explore  the 
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wave  propagation  characteristics  in  the  complex  near-source  environment,  we  conduct 
slowness  analysis  for  these  wavefields.  We  note  that  any  differences  in  wave  excitation 
due  to  model  variations  relative  to  the  reference  flat  surface  case  are  captured  in  the 
differenced  wavefields  and  may  not  be  distinguishable  from  scattering  effects,  but 
statistical  averaging  over  multiple  model  realizations  does  tend  to  isolate  the  scattering 
effects. 

Rg 


Figure  19.  Slowness  analysis  results  in  the  depth-slowness  domain  for  discrete  time 
intervals  as  the  wave  sweeps  through  the  observing  array  for  models  with 
different  source  depths  and  free  surface  parameters.  The  frequency  band  is 
between  1.5  Flz  and  4.5  Flz.  All  the  panels  are  normalized  in  the  same  scale.  In 
each  small  figure,  the  horizontal  coordinate  is  horizontal  slowness  and  the 
vertical  coordinate  is  depth.  The  thick  vertical  line  indicates  the  upper  mantle  S- 
wave  slowness  which  separates  energy  that  leaks  out  of  the  waveguide  (to  the 
left)  from  energy  trapped  in  the  wave  guide  (to  the  right)  that  forms  Lg.  The  PP, 

PS  and  PSM  are  crustal  P-wave,  S-wave  and  upper  mantle  S-wave  slownesses, 
with  their  values  are  0.154,  0.28  and  0.22  s/km,  respectively. 

Slowness  domain  presentation  of  surface  scattering 

To  quantify  the  scattered  wavefield,  we  transform  the  signal  to  the  slowness  domain. 
We  compute  wavefields  for  the  two-layer  velocity  model  listed  in  Table  3.  The  random 
surface  topographic  fluctuation  is  located  above  the  source  and  extends  in  both  directions 
for  20  km.  The  random  topography  has  an  exponential  power  spectrum  and  its  correlation 
distance,  a,  is  0.5  km.  The  rms  fluctuations  used  are  0.15  km  (maximum  peak-to-trough 
0.625  km)  and  0.3  km  (maximum  peak-to-trough  1.281  km).  Either  a  shallow  source 
(depth  =  0.5  km)  or  a  deep  source  (depth  =  3.0  km)  is  used  in  the  simulation.  The 
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synthetic  seismograms  computed  across  a  vertical  array  (see  Figure  17)  are  processed 
using  the  slowness  analysis  method  (for  details,  see  Xie  et  ah,  2005a).  Illustrated  in 
Figure  19  are  energy  distributions  in  a  mixed  horizontal-slowness  and  depth  domain  (see 
Figure  lb),  where  different  rows  are  for  different  combinations  of  source  and  free  surface 
parameters  and  the  successive  frames  in  each  row  are  for  different  time  windows.  The 
time  progression  allows  ready  identification  of  major  phases.  Phases  such  as  Rg,pS,  S* 
and  trapped  energy  from  different  mechanisms  are  labeled  in  the  figure.  Frequency  filters 
can  be  applied  to  the  synthetic  data  before  conducting  the  slowness  analysis.  When  this  is 
viewed  collectively,  we  know  the  energy  distribution  in  combined  domains  of  space, 
slowness,  time  and  frequency.  Hereafter,  we  will  call  this  SPTF  domain  characterization, 
with  the  “F”’  standing  for  slowness. 

Table  4.  Source  and  model  parameters  used  in  numerical  simulations 


Group 

Fixed  parameters 

Variables 

1 

a  =  0.5  km 

Source  depth  =  0.5  km 
Cp  =  infinity 

Qs  =  infinity 

rms  free  surface  fluctuation  = 

0.0,  0.05,  0.10,  0.15,  0.20,  0.25,  0.30,  0.35,  0.40  km 

2 

rms  =  0.15  km 
a  =  0.5  km 

Qy  =  infinity 

Qs  =  infinity 

Source  depth  = 

0.25,  0.50,  0.75,  1.00,  1.50,  2.00,  2.50,  3.00  km 

3 

rms  =  0.15  km 
Source  depth  =  0.5 
km 

Qv  =  infinity 

Qs  =  infinity 

Correlation  length  a  = 

0.4,  0.6,  0.8,  1.0,  2.0,  4.0,  6.0,  8.0,  10.0  km 

4 

a  =  0.5  km 

Source  depth  =  0.5 
km 

gp=  100 

0s  =  50 

rms  free  surface  fluctuation  = 

0.0,  0.05,  0.10,  0.15,  0.20,  0.25,  0.30,  0.35,  0.40  km 

In  Figure  19,  the  solid  vertical  lines  indicate  the  upper  mantle  S'-wave  slowness. 

Wave  energy  to  the  left  of  these  lines  has  incidence  angles  steeper  than  the  critical  angle 
on  the  Moho  and  the  energy  will  leak  to  the  upper  mantle  through  multiple  reflections. 
For  wave  energy  to  the  right  of  these  lines,  total  reflection  will  keep  the  energy  in  the 
crustal  wave  guide,  ultimately  forming  the  Fg-wave  at  long  distances  (e.g.  Frankel,  1989; 
Xie  and  Lay,  1994;  Vogfjord,  1997;  Xie  et  ah,  2005a).  Figure  19a  is  for  a  flat  free  surface 
and  a  relatively  deep  explosion  source.  As  expected,  this  configuration  generates  neither 
noticeable  trapped  energy  nor  clear  i?g-wave.  Figure  19b  is  for  a  flat  free  surface  and  a 
shallow  explosion  source.  We  now  see  the  7?g-wave  developed  at  shallow  depth  and 
trapped  energy  originating  from  the  5* -wave  (Vogfjord,  1997;  Xie  et  al.  2005a).  Figure 
19c  is  for  a  deeper  explosion  source  and  a  free  surface  with  0.15  km  rms  fluctuation. 
Comparing  with  19a,  the  existence  of  a  random  free  surface  generates  significant  trapped 
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energy  from  surfaee  seattering.  Although  the  souree  is  loeated  at  a  depth  of  3.0  km,  Rg 
energy  ean  now  be  seen  at  shallow  depth.  This  enhaneed  i?g-wave  eomes  from  the  free 
surfaee  seattering,  whieh  ean  be  treated  as  shallow  seeondary  sourees.  In  Figures  19d  to 
1 9f,  with  shallower  souree  or  larger  rms  free  surfaee  fluetuations,  a  substantial  amount  of 
trapped  energy  ean  be  generated  from  interaetions  between  the  explosion  souree  and  the 
topographie  fluetuations.  The  energy  partitioned  to  speoifie  regions  in  SPTF  domain  ean 
be  eolleeted  and  used  to  estimate  the  energy  input  into  different  regional  phases.  Using 
different  souree/model  parameters,  the  relationship  between  these  parameters  and  the 
energy  partitioning  proeesses  is  investigated  and  the  underlying  meehanisms  revealed. 
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Figure  20.  Slowness  analysis  results  for  models  with  different  rms  topographic 
fluctuations.  The  source  depth  is  0.5  km.  In  each  time  frame  the  horizontal 
coordinate  is  the  horizontal  slowness  and  vertical  coordinate  is  the  depth.  The  Rg 
energy  is  located  near  the  surface,  with  slowness  similar  to  the  S-wave. 


5.  2.  Contributions  of  Free  Surface  Scattering  to  Cg-Wave  Excitation 

We  adopt  the  three-layer  erust  (Table  2)  as  our  basie  veloeity  model  and  add  random 
free  surface  fluctuations  with  different  statistical  parameters  to  this  basic  model.  The 
random  topography  has  an  exponential  power  spectrum.  It  is  located  above  the  source 
and  extends  in  both  directions  for  20  km.  The  model  geometry  is  similar  to  that  shown  in 
Figure  17,  with  a  vertical  array  composed  of  41x90  receivers  located  between  distances 
30  and  50  km  and  depths  0  and  45  km.  To  isolate  the  effects  of  individual  factors,  we 
vary  individual  parameters  while  keeping  other  parameters  unchanged.  The  varied 
parameters  are  listed  in  Table  4.  Each  model  is  described  by  a  set  of  parameters  including 
the  rms  free-surface  fluctuation,  the  correlation  length  of  the  random  power  spectrum,  the 
source  depth,  and  the  intrinsic  attenuation  (quality  factor  Q).  The  source  depth  is 
measured  locally  by  taking  the  vertical  distance  between  the  source  and  the  rough 
surface.  This  avoids  coupled  variation  of  depth  and  topographic  parameters.  To 
characterize  the  results  statistically,  we  generate  10  realizations  for  each  model.  Synthetic 
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seismograms  are  calculated  for  each  realization  and  processed  separately.  We  then 
average  the  measurements  from  individual  realizations  and  use  their  mean  value  as  the 
final  result  for  a  particular  case. 

5.  3.  Effect  of  free  surface  roughness 

To  investigate  the  effect  of  free  surface  roughness  on  the  energy  partitioning,  we  use 
an  infinite  Q  for  both  P-  and  ^-waves,  a  source  depth  of  0.5  km,  a  correlation  length  of 
0.5  km,  and  vary  the  rms  free  surface  fluctuation  between  0.0  and  0.4  km  (see  group  1  in 
Table  4).  The  slowness  analysis  results  are  illustrated  in  Figure  20  with  different  panels 
for  models  with  different  rms  values.  A  frequency  domain  filter  between  1.5  Hz  and  4.0 
Hz  was  applied  to  the  synthetic  seismograms  before  this  slowness  analysis.  The  Rg 
energy  can  be  seen  at  depths  less  than  3  km  with  a  slowness  similar  to  the  ^-wave.  In 
Figure  20a,  with  rms  =  0,  the  i?g-wave  is  directly  generated  entirely  by  the  explosion 
source.  It  arrives  at  the  receiver  array  between  10  and  12  sec  and  is  labeled  as  “direct  Rg- 
wave”. 

(a)  direct  Rg  energy  (b)  scattered  Rg  energy  (c)  Lg  energy 


.1 


(d)  scattered  loss  of  Rg  energy  (e)  net  scattered  Rg  energy  (f)  scattered  gain  of  Lg  energy 


Figure  21.  Responses  as  functions  of  frequency  and  rms  free  surface  fluctuations.  Top 
row:  the  near  source  responses  of  direct  Rg  ,  scattered  Rg  and  the  Lg-waves,  with 
(a)  ^  (b)  and  (c)  .  Bottom  row:  the  contributions  of  surface 

nRg_direct  _  pRgjiirect  pRg_scatt  _  nRg_scatt 

scattering  to  these  responses,  with  (d)  ^  ,  (e)  ^  and 

(f)  f  .  Note,  a  negative  vertical  coordinate  is  used  in  (d),  where  the  prisms 
with  solid  black  tops  are  below  zero  plane,  while  prisms  with  open  tops  are 
above  the  zero  plane. 

The  presence  of  a  rough  free  surface  causes  scattering  of  different  waves  and 
redistributes  their  energy.  As  shallow  secondary  sources,  the  scattering  generates 
scattered  i?g-wave  which  can  be  observed  in  all  time  windows  in  panels  20b  and  20c.  The 
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shallow  energy  between  8  and  10  see  is  from  body  wave  to  Rg  seattering,  as  it  arrives  at 
the  reeeiver  array  earlier  than  the  direet  Rg.  We  label  this  as  “seattered  i?g-wave”  and  use 
it  to  investigate  body-to-surfaee  wave  seattering.  Due  to  seattering  attenuation,  the  same 
surfaee  fluetuations  that  exeite  the  scattered  Rg  can  also  attenuate  both  direct  and 
scattered  Rg,  generating  scattered  body  waves.  In  panels  20d  and  20e,  due  to  strong 
scattering  from  a  very  rugged  free  surface,  both  direct  and  scattered  waves  are  very 
weak.  Using  panel  20a  as  a  reference,  the  trapped  energy  can  be  quantified  as  a  function 
of  surface  roughness. 

Applying  SPTF  processing  to  the  slowness  analysis,  we  can  separate  the  energy  and 
estimate  the  excitation  of  different  phases  (Xie,  et  ah,  2005a).  A  series  of  band  pass 
filters  is  used  to  obtain  responses  between  1.0  and  4.0  Hz.  Figures  21a  to  21c  illustrate 
the  near-source  response  functions  of  direct  i?g-wave,  scattered  Rg-wave  and  Zg-wave 
(summed  trapped  energy)  as  functions  of  frequency  and  rms  surface  fluctuations.  The  Rg- 
wave  energy  is  obtained  by  picking  energy  peaks  between  depth  0  and  3  km  and 
horizontal  slowness  0.18  -  0.40  s/km.  The  Tg-wave  energy  is  picked  from  peaks  between 
depths  3  and  45  km,  and  slowness  0.22  -  0.40  s/km.  The  vertical  coordinate  is  the  square 
root  of  normalized  energy.  The  source  spectrum  has  been  removed  from  the  normalized 
energy  and  the  results  are  equivalent  to  the  response  functions  in  equations  (11)  and  (12). 
Since  ii  oc  ,  where  A  is  the  displacement  amplitude,  the  square  root  of  energy  can  be 
compared  to  the  wave  amplitude  after  being  scaled  with  /'' .  For  the  response 
function  ,  the  energy  is  mainly  located  at  low  frequencies  and  drops  quickly  with 

increasing  rms  values.  Note  that  the  histogram  labeled  with  rms  =  0  indicates  the 
response  for  the  flat  earth  model.  The  net  contribution  of  the  surface  scattering  to  direct 
i?g-wave  can  be  obtained  by  subtracting  the  fiat  earth  response  from  the  total  response 
(Figure  2 Id).  The  generally  negative  values  indicate  the  energy  loss  from  direct  Rg  due  to 
the  scattering  attenuation,  which  is  proportional  to  the  rms  values.  Note  that  a  negative 
vertical  scale  is  used  in  Figure  2 Id,  and  the  prisms  with  solid  black  tops  are  plotting 
below  the  plane  zero,  while  prisms  with  open  tops  are  plotting  above  the  plane.  At  very 
low  rms  values,  the  i?g-wave  gains  some  energy,  which  appears  to  be  caused  by  body-to- 
surface  wave  scattering  providing  slightly  more  Rg  energy  than  the  energy  loss  from  Rg 
due  to  scattering. 

The  time  windows  of  direct  and  scattered  Rg-waves  partially  overlap.  The 
application  of  a  narrow-band  filter  causes  broadening  of  the  Rg  impulse  and  results  in 
some  leakage  of  direct  Rg-wave  energy  into  the  scattered  Rg  time  window  (see  first  row 
in  Figure  21b).  By  subtracting  this  energy,  we  obtain  the  net  scattering  contribution, 
j^Rg  scatt  _  ^  jg  shown  in  Figure  21e.  This  energy  increases  with  increasing 

rms  up  to  a  moderate  rms  value,  then  decreases  with  further  increasing  rms.  This  suggests 
that  small  to  moderate  topographic  fluctuations  can  excite  secondary  i?g-wave,  but  a  very 
bumpy  free  surface  will  dramatically  attenuate  the  short  period  i?g-wave  as  a 
consequence  of  both  weakened  excitation  and  strong  scattering.  Both  low  and  high 
frequency  scattered  i?g-waves  can  be  identified.  The  high  frequency  content  falls  faster 
with  increasing  rms  values  than  its  lower  frequency  counterpart. 
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Shown  in  Figures  21c  and  21f  are  the  Zg-wave  response,  (i.e.,  summed  trapped 
energy  in  the  waveguide),  and  net  scattering  contribution  to  Zg-wave,  -Rjf , 
respectively,  as  functions  of  surface  topography  rms  value.  Comparing  Figures  21d  and 
2 If,  at  low  frequencies,  the  maximum  scattered  loss,  ,  from  Rg  is  0.04,  while  the 

maximum  scattered  gain  of  Lg  is  0.16  at  rms  value  0.40  km.  Assuming  that  all  scattered 
loss  of  Rg  goes  to  Lg,  there  still  needs  to  be  some  Lg  energy  contributed  from  other 
sources,  e.g.,  body  wave  scattering  at  the  rough  free  surface.  A  broadband  C-wave  can 
excite  low  frequency  S'-waves  because  the  surface  scattering  serves  as  shallow  secondary 
sources.  Our  method  cannot  completely  disentangle  the  percentage  of  scattered  Rg  going 
into  Lg,  but  the  evidence  for  significant  low  frequency  contribution  is  (1)  the  wavefield 
snapshot  (Figure  18d)  shows  considerable  scattered  Rg  going  into  trapped  waves,  (2)  the 
spectrum  shape  of  energy  loss  of  Rg  is  similar  to  the  energy  gain  of  Lg  (Figures  21  d  and 
2 If),  and  Lg  energy  gain  shows  strong  source  depth  dependence,  as  shown  below.  At 
high  frequencies,  Rg  is  weakly  excited,  so  its  scattering  provides  little  energy  for  Lg- 
waves:  most  of  the  high  frequency  scattered  Lg  energy  comes  from  body  wave  scattering. 


scattered  Rg-wave 


(a)  local  source  depth=0.25  km 


direct  Rg-wave 


S-s 

T35 


^  1 

f 

‘  ■  i 

9  NT 

10 

; 

10  ft 

"ir 

12.ft  •  1 

13  06  j  • 

'1'- 

i 

i6«j: 

< 

J 

1  • 

Ua 

J;  i 

'_4 

- 

'  1  .  1 
«:  ! 

1  8.01  '  :  '  90* 

1 

1  ! 

10^ 

: 

'  ' 

V 

■  trapped  energy 


'X^ 

13  05 

14  Os !  : 

•  1  • 

15.Cs 

; 

A 

4 

:  : 

:  1  ■ 

’  J*. 

m 

(c)  local  source  depth=1  Okm 


£8 


(d)  local  source  depth=2.0km 


£8 

&8 


(e)  local  source  depth=3.0km 


£8- 

g-S 

•D@J 


'T^l  '  :  '  90^ 

i  1  :  1 

1  95|'|  10* 

lO.ft  1  ■  ~ 
1  1 

12  ft 

* 

13ft  i  : 

:  i ; 

14  ft 

15%  I  ; 

hli .. 

1  '  1  ;  j. 

1  A;  i;  1 

>  j  ^ 

t 

1, 

4 

1  m 

: 

:  i  : 

:  i  i 

ec&j  : 

’ '  9  06  j  : 

1 9  T  i : 

10  ft 

■ 

12ft  1  t 

13ft  i  :  1 '  14ft  {  ; 

15ft 

I  : 

I  > 

i 

‘  f,  ; 

sT 

r 

M 

!*|  : 

••  1  ; 

:  1  ; 

T 

;  1  ;  1  :  : 

.1:  '■  • 

• 

i  1 

T  60%  ■  904  ’  ■ 

:  '  :  1  .  i  ; 

9.5f  1  '  1  10  ft  1 

♦  .  1  " 

10  ft  1  L 

11  ft  '  ^ 

12ft  1  : 

lOOsj  :  ,  t4ft  1  :  1 

15  ft  1  :  1 

~  6 : 1 

1  ■!; 

s 

Hi  ! 

-* _ .  ft:  1  . 

■  L:  :  .  • .  1 

■  1 

Figure  22.  Similar  to  Figure  20  except  each  panel  is  for  a  different  source  depth. 

5.  4.  The  Effect  of  Source  Depth 

To  investigate  the  effect  of  source  depth  on  the  energy  partitioning,  we  fix  the  rms 
free  surface  fluctuation  at  0.15  km,  the  correlation  length  as  0.5  km  and  vary  the  source 
depth  between  0.25  km  and  3.0  km  (see  group  2  of  Table  4).  The  slowness  analysis 
results  are  shown  in  Figure  22  with  each  panel  being  for  a  model  with  a  particular  source 
depth.  A  prominent  feature  is  that  shallower  sources  generate  more  Rg  and  Lg  energy. 
Applying  the  SPTF  processing  to  the  slowness  analysis,  we  obtain  energy  response 
spectra  for  different  phases.  Figures  23a  to  23c  present  the  response  functions  of  direct 
Rg- ,  scattered  Rg-,  and  Tg-waves  with  respect  to  frequencies  and  source  depths.  With 
increasing  source  depth,  the  Rg  energy  falls  quickly.  Shown  in  Figures  23  d  to  23  f  are 


29 


separated  net  contributions  from  the  free  surface  scattering  obtained  by  subtracting  the 
energy  of  the  flat  earth  model  from  that  for  a  random  free  surface.  For  direct  Rg  from  a 
shallow  source,  a  large  amount  of  energy  is  lost  due  to  scattering.  However,  for  deeper 
sources,  the  scattering  adds  a  small  amount  of  energy  to  the  i?g-wave  (Figure  23d).  The 
scattered  Lg-wave  shows  an  apparent  increase  for  sources  shallower  than  about  0.5  km 
which  suggests  a  surface  wave  origin.  The  energy  budget  between  Figures  23d  and  23f 
are  also  comparable  for  lower  frequencies  and  shallow  sources.  For  Tg-wave  with 
frequencies  of  3  Hz  or  higher,  the  energy  should  come  from  body  wave  surface 
scattering. 

For  deeper  sources,  the  lower  frequency  Tg-wave  is  very  weak,  indicating  that  both 
Rg-to-Lg  and  body-to-Tg  conversions  are  weak  within  this  frequency-depth  range. 


(a)  direct  Rg  energy 


(b)  scattered  Rg  energy 


(c)  Lg  energy 


Figure  23.  Similar  to  Figure  21,  except  responses  are  functions  of  frequencies  and  source 
depths.  The  rms  surface  fluctuation  is  0.15  km  ,  with  correlation  length  0.5  km. 


5.  5.  The  Effect  of  Correlation  Length 

To  investigate  the  effect  of  correlation  length  on  the  energy  partitioning,  we  fix  the 
rms  free  surface  fluctuation  as  0.15  km,  the  source  depth  at  0.5  km  and  vary  the 
correlation  length  between  0.4  km  and  10  km.  These  parameters  are  listed  in  group  3  of 
Table  4.  The  slowness  analysis  results  are  shown  in  Figure  24  with  each  panel  for  a 
model  with  different  correlation  length.  The  response  functions  of  Rg-,  scattered  Rg-  and 
Tg-waves  are  shown  in  Figures  25  a  to  25  c.  The  horizontal  coordinate  is  frequency  and 
different  rows  are  for  different  correlation  lengths.  From  these  response  functions,  we  see 
that  the  last  row  (for  a  correlation  length  of  10  km)  is  almost  the  same  as  that  for  a  flat 
earth  model.  This  indicates  that  a  very  smooth  long  wavelength  free  surface  fluctuation 
has  almost  no  effect  on  these  waves.  Figures  25 d  to  25 e  are  net  energy  loss  or  gain  from 
scattering.  For  models  with  correlation  length  shorter  than  4  km,  the  surface  scattering 
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apparently  contributes  to  the  generation  of  the  trapped  energy  in  the  waveguide  (Figure 
25 f).  For  very  long  correlation  lengths,  the  random  free  surface  behaves  more  like  a  flat 
free  surface. 


scattered  Rg-wave 


direct  Rg-wave 
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Figure  24.  Similar  to  Figure  4  except  each  panel  is  for  a  model  with  different  correlation 
distance. 


In  Figure  25  f,  the  net  scattered  contribution  to  the  Zg-wave,  the  response  function 
falls  with  increase  of  the  correlation  length  at  all  frequencies.  For  frequencies  used  in  the 
simulation  (1  to  4  Hz)  and  ^-wave  velocity  in  the  top  layer  (3.2  km/sec),  the  wavelengths 
are  between  0.8  and  3.2  km.  We  calculate  the  response  as  a  function  of  normalized  scale 
length,  ka  ,  where  k  =  In/X ,  2  is  the  wave  length  and  a  is  the  correlation  length,  and 
present  the  behavior  in  Figure  26.  The  maximum  scattering  happens  around  ka  =  1  and 
decreases  for  larger  .  Extension  of  the  calculation  to  smaller  ka  is  limited  by  the  grid 
size  used  in  the  boundary  element  calculation  and  the  dimension  of  the  receiver  array. 

5.  6.  The  Effect  of  Intrinsic  Attenuation 


Scattering  from  topographic  fluctuations  occurs  in  the  uppermost  crust,  which  is 
usually  a  low  Q  layer.  In  addition,  the  scattering  increases  the  propagation  distances, 
especially  for  high  frequency  waves.  Attenuation  will  thus  strongly  affect  the  scattering 
and  the  overall  energy  distribution  of  an  explosion  source.  To  assess  the  effects  of 
shallow  attenuation,  we  use  a  set  of  model  parameters  similar  to  that  used  for  testing  the 
effect  of  rms  fluctuations  except  we  replace  the  infinite  Q  in  the  top  10  km  with  Qp=  100 
and  Qs  =  50.  These  parameters  are  listed  in  group  4  of  Table  4  and  the  results  are  shown 
in  Figure  27.  Comparing  Figure  27  with  Figure  21,  two  prominent  features  can  be 
identified.  First,  compared  to  the  purely  elastic  case,  there  is  significant  energy  loss  in  the 
model  with  intrinsic  attenuation.  For  example,  the  maximum  amplitude  (square  root 
energy)  drops  approximately  35%  for  direct  Rg-wavc,  35%  for  scattered  7?g-wave  and 
40%  for  Tg-wave.  Second,  the  short  period  waves  experience  even  larger  attenuation  than 
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long  period  waves.  This  is  espeeially  true  for  the  seattered  Rg-wave  and  the  Tg-wave.  By 
using  low  Q  values  in  the  ealeulation,  our  results  should  give  fairly  extreme 
eharaeterization  of  the  effeet  of  attenuation. 


(f)  scattered  gain  of  Lg  energy 


(a)  direct  Rg 


(b)  scattered 


(d)  scattered  loss  of  Rg  energy 


(e)  net  scattered 


Figure  25.  Similar  to  Figure  21,  except  responses  are  functions  of  frequencies  and 

correlation  lengths. 

5.  7.  Variance  of  Statistical  Results 

To  obtain  statistical  relationships  between  the  source/model  parameters  and  the 
energy  partitioning  proeess,  we  average  the  measurements  from  individual  realizations. 
To  investigate  the  varianee  from  a  group  of  realizations  with  the  same  statistieal 
parameters,  we  compare  results  for  models  with  rms  value  0.15  km,  eorrelation  length  0.5 
km,  souree  depth  0.5  km  and  infinite  Q.  The  response  funetions  for  different  phases  are 
presented  in  Figures  28a  to  28o,  where  the  horizontal  coordinates  are  frequeney  and 
different  rows  are  for  different  realizations.  The  response  eurves  for  different  realizations, 
their  mean  values  and  standard  deviations  are  also  shown  in  Figures  28d  to  28f.  The 
primary  Rg-wave  is  relatively  stable.  The  scattered  Rg-wave  has  large  variations  at 
frequencies  below  about  2  Hz,  while  the  Tg-wave  shows  large  variations  at  frequencies 
above  2  Hz.  Although  based  on  the  same  statistieal  parameters,  response  functions  from 
different  realizations  show  variations  in  amplitudes  and  loeal  minima,  suggesting  that  the 
partitioning  is  partially  affeeted  by  deterministie  features  very  elose  to  the  source. 
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ka 

Figure  26.  Net  scattered  Lg-energy  as  a  function  of  normalized  scale  factor  ka. 


(a)  direct  Rg  energy 


(b)  scattered  Rg  energy 


(c)  Lg  energy 


Figure  27.  Similar  to  that  shown  in  Figure  21,  except  a  low  Q  top  layer  is  used  in  the 
simulation.  To  facilitate  comparison,  the  same  vertical  scale  as  in  Figure  21  is 
used  here. 
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(a)  direct  Rg  energy 


(b)  scattered  Rg  er>ergy 


(c)  Lg  energy 


Figure  28.  Top  row:  Responses  as  funetions  of  frequencies  and  different  random 
realizations,  for  (a)  direct  Rg-wave,  (b)  scattered  Rg-wave  and  (c)  Lg-wave.  The 
bottom  row  shows  response  spectra  from  individual  realizations,  their  mean 
values  and  standard  deviations. 


Exptosion  30  km 


Figure  29.  Cartoon  showing  the  configuration  of  the  3D  velocity  model,  source,  and 
receiver  array. 

6.  REGIONAL  PHASE  EXCITATION:  SCATTERING  IN  THREE- 
DIMENSIONS 


6. 1.  The  3D  Volumetric  Scattering  and  Excitation  of  SH  Component 

Although  most  Lg  observations  are  taken  from  the  vertical  component  of 
seismograms,  the  tangential  component  ohen  has  as  much  energy  as  the  vertical 
component  (e.g.,  Stevens,  et  al.  2003).  Some  researchers  pointed  out  that  at  close  range, 
clear  SH  energy  comparable  to  that  on  the  SV  component  can  be  observed  and  must  be 
generated  in  the  source  region.  Any  Lg-wave  excitation  theory  must  provide  an 
explanation  for  these  observations.  Since  2D  geometry  decouples  the  P-SV  problem  and 
the  SH  problem,  it  does  not  provide  any  information  on  the  coupling  between  the  source 
and  SH  component.  We  have  developed  the  3D  local  slowness  analysis  method.  A  three 
dimensional  velocity  model  is  used  to  simulate  the  near  source  environment.  The  size  of 
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the  model  is  30  kmx  30  kmx  20  km,  and  the  upper  erust  structure  from  the  EK  model  is 
used  as  the  background  velocity.  To  test  the  effect  of  the  heterogeneities  on  the  P-S 
coupling,  7%  RMS  broadband  random  perturbation  is  added  to  P-  and  S-wave  velocities 
and  the  density  between  two  cylindrical  surfaces  around  the  z-axis.  An  eighth-order, 
staggered  format  3D  elastic  finite  difference  code  is  used  to  generate  synthetic 
seismograms.  The  grid  interval  used  is  0.1  km.  The  explosion  source  is  located  at  depth 
0.5  km  and  the  dominant  frequency  is  about  3  Hz.  A  10x10x30  3D  receiver  array  is 
located  at  epicentral  distance  28  km  and  azimuth  direction  45  degree  (see  Figure  29  for 
model  configuration). 

Distance  km  Distance  km 

0  5  10  15  20  25  30  0  5  10  15  20  25  30 


Firgure  30.  Wavefield  snapshot  for  the  layered  background  model  (left)  and  a  laterally 
heterogeneous  model  with  7%  RMS  random  fluctuations.  Shown  here  is  the 
vertical  component  of  the  displacement  field. 
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Figure  31.  Comparison  between  synthetic  seismograms  and  energy  distribution  in 
horizontal  slowness  domain  for  the  layered  velocity  model  (left  panel)  and 
velocity  model  with  7%  RMS  velocity  fluctuations  (right  panel).  Receivers  are 
located  at  depth  1  km.  All  three  components  of  the  seismograms  are  normalized 
jointly  but  the  slowness  distribution  for  the  radial  component  of  the  P-wave  has 
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been  multiplied  by  a  factor  of  0.1.  The  P  and  Rayleigh  waves  can  be  clearly  seen 
from  the  radial  and  vertical  components.  Note  that  the  tangential  component  for 
the  layered  model  is  zero.  For  the  random  velocity  model,  scattered  waves  in 
both  the  tangential  and  vertical  components  can  be  seen.  Note  that  much 
scattered  energy  falls  outside  of  the  upper  mantle  slowness  (dash  line  circle)  and 
can  be  trapped  into  the  crustal  wave  guide  to  form  the  Lg  wave. 

To  investigate  the  effeet  of  the  random  velocity  perturbations,  models  both  with  and 
without  the  random  fluctuations  are  calculated  and  their  wavefield  snapshots  are  shown 
in  Figure  30.  No  surface  roughness  is  included  here.  Scattered  energy  appears  in  the 
model  with  random  velocity  fluctuations.  Shown  in  Figures  31  are  synthetic  seismograms 
and  the  2D  horizontal  slowness  analysis  from  selected  time  windows  for  both  the 
background  and  random  models.  The  source  depth  is  0.5  km  and  receiver  arrays  are 
located  at  depth  1  km.  For  the  background  model  (left  panel),  the  P  and  Rayleigh  waves 
can  be  clearly  seen.  Note  the  tangential  component  is  zero.  All  three  components  of  the 
seismograms  are  normalized  jointly,  but  the  slowness  distribution  for  the  radial 
component  of  the  P-wave  has  been  multiplied  by  a  factor  of  0. 1 .  For  the  random  model 
with  7%  RMS  velocity  fluctuations  (right  panel),  there  are  scattered  waves  in  all  three 
components  including  the  tangential  (SH)  component.  Much  of  the  scattered  energy  falls 
outside  of  the  upper  mantle  slowness  (dashed  line  circle)  and  can  thus  be  trapped  into  the 
crustal  waveguide  to  form  the  regional  phases. 


Figure  32.  The  wave  energy  distribution  in  the  2D  horizontal  slowness  domain  for  the 
vertieal  (left)  and  tangential  (right)  eomponents.  The  receiver  array  is  located  at 
epicentral  distance  28  km  and  depth  7  km.  Note  the  energy  distribution  outside  of 
the  upper  mantle  S-wave  velocity. 

Figure  32  gives  the  2D  horizontal  slowness  analysis  for  vertical  and  tangential 
components  for  the  random  velocity  model.  The  receiver  array  is  located  at  depth  7  km 
and  epicentral  distance  28  km.  In  the  left  panel  for  the  vertical  component,  we  can  see 
that  both  P-wave  (around  5.2  s  to  5.6  s)  and  pS-wave  (around  6.8  to  7.2  s)  share  the  same 
slowness  and  will  not  contribute  to  the  guided  wave  mode.  Within  the  time  window 
between  7.2  and  9.2  s,  there  is  scattered  P-SV  energy  propagating  with  a  larger  apparent 
horizontal  slowness.  The  energy  falls  beyond  the  upper  mantle  S-wave  slowness  and  will 
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be  trapped  in  the  waveguide  to  form  the  P-SV  Lg  wave.  In  the  right  panel  for  the 
tangential  component,  there  is  considerable  SH  type  scattered  energy  between  7.2  and  9.2 
s  with  slowness  greater  than  the  upper  mantle  S-wave  slowness.  The  following  features 
can  be  found  from  the  results:  (1)  There  is  considerable  SH  energy  excited  due  to  the  near 
source  lateral  heterogeneities;  (2)  The  SH  energy  appears  to  be  generated  through  P-pS- 
SH  or  SV-SH  since  the  SH  component  is  relatively  weak  immediately  following  the  P- 
wave  front.  The  scattered  energy  has  a  broader  azimuth  range  than  the  direct  arrivals, 
which  may  provide  a  hint  to  reveal  the  actual  scattering  mechanisms.  (3)  Both  SV  and  SH 
scattered  energy  can  fall  into  the  proper  slowness  region  and  form  the  crustal  guided 
wave. 


Figure  33.  The  3D  velocity  model  including  both  irregular  free  surface  and  volumetric 
heterogeneity  (shown  here  is  one  quarter  of  the  model).  The  rough  topography  is 
located  between  distances  5  km  and  15  km.  The  volumetric  heterogeneity  is 
located  between  distances  5  km  and  1 5  km,  and  between  depths  1  km  and  5  km. 
The  synthetic  seismograms  are  collected  from  a  3D  array  located  at  epicentral 
distance  28  km. 


6.  2.  The  3D  Simulation  Including  both  Rough  Free  Surface  Scattering 
and  Volumetric  Scattering 

To  include  both  volumetric  heterogeneity  and  irregular  topography  in  a  3D  simulation 
and  investigate  their  effects  on  near  source  energy  partitioning,  we  use  the  3D  finite- 
difference  method  developed  by  Zhang  and  Chen  (2006).  In  this  finite-difference 
algorithm,  the  physical  domain  is  discretized  by  boundary-conforming  grids.  The 
irregular  surface  is  transformed  into  a  ‘fiat’  surface  in  computational  space,  thus  avoiding 
the  artifact  usually  caused  by  using  the  staircase  approximation.  To  satisfy  the  free 
surface  boundary  conditions,  the  method  uses  the  stress  image  method  on  the  irregular 
topography.  With  this  method,  approximately  10  grid  points  per  shortest  wavelength  is 
usually  enough  to  maintain  the  global  accuracy  of  the  simulation.  For  simulation,  we  use 
a  60km  X  60km  x  30km  velocity  model  (grid  size  600  x  600  x  300  ).  Shown  in  Figure  33  is  one 
quarter  of  the  model.  An  explosion  source  is  located  at  the  center  of  the  model  with  a 
depth  of  0.5  km.  The  rough  topography  is  located  between  epicenter  distances  5  and  15 
km.  The  volumetric  heterogeneity  is  located  between  horizontal  distances  5  and  15  km 
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and  between  depths  1  km  and  5  km.  A  xiv^,  =10x10x16  ( 4kmx4kmx30km) 

receiver  array  located  at  epicenter  distance  28  km  is  used  to  collect  the  synthetic 
seismograms. 


x-component 


y-component 


z-component 


(d) 

Layer 
+  volumetric 
+  topography 


Figure  34.  Snap  shots  for  an  explosion  souree  in  different  veloeity  models.  Rows  (a)  to 
(d)  are  for  a  layered  model,  a  layered  model  with  volumetrie  random  veloeity 
fluetuation,  a  layered  model  with  irregular  topography,  and  a  layered  model  with 
both  volumetrie  heterogeneity  and  random  topography.  The  three  eolumns  are  for 
X-,  y-  and  z-eomponents. 


Shown  in  Figure  34  are  snapshots  for  an  explosion  source  in  different  velocity 
models.  From  rows  (a)  to  (d)  are  for  a  layered  model,  a  layered  model  with  volumetric 
random  velocity  fluctuation,  a  layered  model  with  irregular  topography,  and  a  layered 
model  with  both  volumetric  heterogeneity  and  random  topography.  The  three  columns 
are  for  x-,  y-  and  z-components.  The  rough  surface  has  an  exponential  power  spectrum,  a 
rms  fluctuation  of  0.2  km,  and  a  correlation  distance  of  2.0  km.  The  volumetric 
heterogeneity  has  an  exponential  power  spectrum,  a  rms  perturbation  of  5%,  and  a 
correlation  length  of  2.0  km.  For  a  layered  model  shown  in  row  (a),  an  explosion  source 
generates  P,  pS,  Rg  and  some  interface  reflections,  which  are  all  P-SV  waves.  As 
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examples,  we  see  the  x-eomponent  is  zero  in  plane  jOz  and  y-eomponent  is  zero  in  plane 
xOz.  For  models  with  volumetrie  and/or  rough  topography,  i.e.,  rows  (b)  to  (d),  there  are 
apparent  seattered  waves  following  the  primary  phases.  Partieularly,  SH-waves  appear  in 
tangential  eomponent,  e.g.,  the  non-zero  x-eomponent  inyOz  plane  and  non-zero  y- 
eomponent  in  xOz  plane. 


Shown  in  Figure  35  are  synthetie  seismograms  and  slowness  analysis  for  (a)  a 
layered  erustal  model  and  (b)  a  layered  model  plus  free  surfaee  and  volumetrie 
heterogeneity.  The  model  geometry  is  shown  in  Figure  33.  The  top,  middle  and  bottom 
rows  are  for  radius,  tangential  and  vertieal  eomponents,  respeetively.  The  rough  surfaee 
has  an  exponential  power  speetrum,  a  rms  fluetuation  of  0.2  km  and  a  eorrelation  distanee 
of  2.0  km.  The  volume  heterogeneity  has  a  5%  veloeity  perturbation,  an  exponential 
power  speetrum  and  a  eorrelation  length  of  2  km.  The  synthetie  seismograms  are  filtered 
between  0.5  and  5.0  Hz  and  normalized  using  the  same  seale.  The  horizontal  distanee 
from  the  souree  to  the  eenter  of  the  reeeiver  array  is  about  28  km.  The  2D  slowness 
analyses  are  eondueted  in  horizontal  plane  (the  map  view)  for  phases  P,  PP,  PS,  Rg  and 
eoda  waves.  The  results  show  that  for  an  explosion  souree  in  a  layered  model,  there  is  no 
SH  wave  and  the  energy  is  distributed  within  the  upper  mantle  S  wave  slowness  (the 
dashed  eirele).  For  a  model  with  both  rough  free  surfaee  and  the  volumetrie 
heterogeneity,  the  primary  P  and  Rg  waves  are  strongly  attenuated  due  to  the  seattering. 
From  slowness  analysis,  there  are  both  SV  and  SH  waves  with  their  slowness  distributed 
outside  the  upper  mantle  S-wave  slowness,  indieating  part  of  their  energy  has  been 
trapped  in  the  erustal  waveguide  to  form  the  Lg-wave. 

6.3.  The  Lg  excitation  Caused  by  the  3D  Surface  Scattering 

We  further  investigate  the  effeet  of  3D  surfaee  seattering  on  the  exeitation  of  Lg 
wave.  We  use  a  model  eonfiguration  similar  to  that  shown  in  Figure  33.  The  eorrelation 
length  of  the  surfaee  random  topography  is  2.0  km  and  the  rms  fluetuation  is  0.15  km.  A 
2.5  Hz  Rieker  wavelet  is  used  in  the  simulation.  We  first  eonduet  the  loeal  slowness 
analysis  and  ealeulate  the  energy  distribution  in  horizontal  slowness  domain  at 

different  depths  (similar  to  that  shown  in  Figure  Id).  Then  we  integrate  the  energy  in  the 
(p:„Py)  domain  along  the  equal  slowness  eireles  and  obtain  the  energy  distribution  as  a 

funetion  of  horizontal  slowness  and  depth.  We  ealeulate  and  eompare  the  energy 
distribution  for  the  layered  EK  model,  EK  model  plus  Gaussian  random  topography  and 
EK  model  plus  exponential  random  topography. 

Shown  in  Eigure  36  is  for  the  EK  model.  Similar  to  that  in  the  2D  ease  (see 
Seetion  5),  in  a  model  with  high  veloeity  top  layer,  the  energy  of  pP  wave  is  mostly 
distributed  to  the  left  of  the  upper  mantle  S-slowness.  The  weak  Eg  wave  (trapped 
energy)  is  generated  from  the  S*  phase.  The  Rg  wave  ean  be  seen  at  the  shallow  part  of 
the  model.  There  is  no  energy  partitioned  to  the  tangential  eomponent.  Shown  in  Eigure 
37  is  for  the  EK  model  plus  the  Gaussian  random  topography.  We  see  part  of  the  energy 
is  transferred  from  P-SV  mode  to  the  SH  mode.  A  Gaussian  random  surfaee  model  is 
mostly  eomposed  of  smoothly  varying  fiuetuations  whieh  generate  weak  seattering.  The 
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Rg  wave  is  weakened  but  the  Love  type  wave  is  generated  at  the  shallow  crust.  Shown  in 
Figure  39  is  for  the  EK  model  plus  the  exponential  random  topography.  Comparing  to  a 
Gaussian  random  model,  an  exponential  random  model  contains  more  small  scale 
topographic  fluctuations  which  tends  to  deflect  more  energy  towards  wide-angle 
scattering,  thus  is  more  effective  in  transmit  energy  in  the  slowness  domain.  The 
scattered  energy  can  be  clearly  observed  following  the  pS  wave  front  in  both  P-SV  and 
SH  mode.  There  is  more  energy  being  shifted  across  the  upper  mantle  S-slowness  and 
being  trapped  in  the  waveguide.  Generally,  a  model  with  more  small  scale  topographic 
fluctuations  tends  to  generate  more  trapped  energy. 


time  (s)  time  (s) 


(a)  (b) 

Figure  35.  Synthetic  seismograms  and  slowness  analysis  for  (a)  a  layered  crustal  model 
and  (b)  a  layered  model  with  rough  free  surface  and  volumetric  heterogeneity. 

The  model  geometry  is  shown  in  Figure  33. 

By  applying  time  windows  to  the  successive  time  frames  in  Figures  36-38  and 
summing  up  the  energy  within  time  windows,  we  obtain  energy  for  different  phases  as 
well  as  for  different  components.  Table  5  lists  the  normalized  energy  for  Rg  and  Eg 
waves  from  different  topographic  models.  We  see  that  the  layered  EK  model  has  large  Rg 
energy  but  very  weak  Eg  energy.  The  energy  for  tangential  component  is  zero.  For 
models  with  random  topographic  fluctuations,  the  Rg  wave  loses  energy  but  the  Eg  wave 
gains  energy.  There  is  energy  transferred  from  the  P-SV  mode  to  the  SH  mode  as  well. 

An  exponential  random  model,  which  has  more  small  scale  fluctuations,  appears  to 
provide  more  efficient  coupling  compared  to  a  Gaussian  random  model.  The  above 
results  may  suggest  that  part  of  the  Rg  energy  is  transferred  into  the  Eg  energy  through 
the  surface  scattering.  However,  the  body  wave  energy,  which  has  not  been  investigated 
here,  may  also  contribute  to  this  process. 
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depth  depth  depth  depth  depth 


(c)Tangential  component 


EK  +  Gaussian  topography 


Figure  37.  Wave  energy  in  horizontal  slowness  and  depth  domain  in  EK  model  plus  free 
surface  rough  topography.  The  random  free  surface  has  a  Gaussian  power 
spectrum. 
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Figure  38.  Wave  energy  in  horizontal  slowness  and  depth  domain  in  EK  model  plus  free 
surface  rough  topography.  The  random  free  surface  has  an  exponential  power 
spectrum. 


Table  5.  Normalized  energy  E/Eq  for  Rg  and  Lg  waves  in  different  random  surfaee 
models. 


Waves  and  components 

EK  model 

EK  model  plus  Gaussian 
topography 

EK  model  plus 
exponential  topography 

H 

5.4x10^'^ 

1.1x10^® 

1.6x10^'^ 

Rg 

V 

2.3x10^*^ 

4.0x10^'’ 

7.0x10^'’ 

T 

0 

2.8  xlO^** 

1.0x10^*^ 

Energy 

2.8x10^* 

Compare  to  EK 

-2.0x10^*^ 

-9.0x10^“^ 

H 

2.8x10^'’ 

5.5x10^'’ 

4.3x10^^ 

Lg 

V 

3.5x10^'^ 

1.1x10^* 

8.6  xlO^'^ 

T 

0 

1.1x10^® 

2.0x10^^ 

Energy 

6.3x10^® 

Compare  to  EK 

+1.1x10^*^ 

+2.7x10^^ 

42 


7.  CONCLUSION 


We  use  2D  and  3D  finite-differenee  method  and  2D  boundary-element  method  to 
simulate  the  seattering  in  the  near-souree  environment.  A  loeal  slowness  analysis  method 
is  developed  to  investigate  near-souree  energy  partitioning  and  Zg-wave  exeitation  for 
explosive  sourees.  The  method  has  two  major  advantages.  First,  it  allows  us  to  study  the 
near-souree  proeesses  in  multiple  domains  including  space,  time,  slowness  and 
frequency.  This  provides  an  opportunity  to  isolate  different  mechanisms  within  the 
complex  near-souree  environment.  Second,  this  method  can  be  applied  at  a  close  range, 
well  before  the  Tg-wave  is  actually  formed.  It  provides  us  with  uncontaminated  near- 
souree  information  by  calculating  a  relatively  small  velocity  model  with  very  fine  near- 
souree  structures.  Since  this  is  a  very  efficient  method,  we  can  use  it  to  investigate  a 
broad  frequency  band  and  to  test  a  large  number  of  source-model  parameters. 

We  investigate  the  contributions  of  P-pS-io-Lg  conversion  and  5'*-to-Tg 
excitation  using  2D  finite-difference  method  and  models  with  near-souree  random 
velocity  fluctuations.  The  results  show  that  the  near  source  scattering  plays  an  important 
role  in  generating  the  Lg  wave  from  an  explosion  source.  The  contribution  of  S*-io-Lg  is 
concentrated  at  low  frequencies  and  occurs  for  very  shallow  source  depths.  The 
contribution  of  P-pS-io-Lg  coupling  in  the  presence  of  near-souree  small-scale  random 
heterogeneities  is  concentrated  at  high-frequencies.  The  excitation  spectra  of  these 
mechanisms  are  calculated. 

We  use  the  2D  boundary-element  simulation  to  investigate  the  effect  of 
topographic  scattering  on  the  near-souree  energy  partitioning  for  an  explosion  source. 
Random  topographic  models  with  different  statistical  properties,  variable  source  depth, 
and  different  Q  models  are  included  in  the  investigation.  The  responses  of  different 
phases  as  functions  of  frequency  and  source/model  parameters  are  calculated  and  their 
energy  budget  evaluated.  The  results  reveal  that  free  surface  scattering  has  strong  effect 
on  near  source  energy  partitioning.  The  scattering  process  can  excite  the  Rg-wave  for  a 
moderately  rugged  topography,  but  dramatically  attenuates  short  period  Rg-waves  when 
the  surface  becomes  too  rugged.  For  models  with  a  high  velocity  shallow  crust,  the  free 
surface  scattering  provides  an  important  mechanism  that  transfers  energy  for  an 
explosion  source  into  the  Lg-wave  in  the  near  source  region.  At  lower  frequencies  and  for 
a  moderately  rugged  free  surface,  the  Rg-to-Lg  transfer  is  relatively  efficient.  At  higher 
frequencies  and  for  a  very  rugged  free  surface,  the  body-to-Lg  transfer  may  dominate  the 
process.  The  correlation  length  of  the  random  free  surface  fluctuation  provides  specific 
frequency  dependence  to  the  transfer  function,  with  maximum  coupling  near  ka=\. 
Intrinsic  attenuation  within  the  uppermost  crust  has  a  strong  effect  on  the  energy  transfer 
through  surface  scattering,  with  high  frequency  content  losing  energy  faster  than  the 
lower  frequency  waves. 

We  expand  our  simulation  and  analysis  to  fully  3D  models  with  both  volumetric 
heterogeneity  and  rough  free  surface.  Investigating  the  near  source  process  in  3D  is  more 
realistic  because  it  gives  correct  scattering  pattern  and  geometric  spreading.  In  addition,  a 
general  3D  model  provides  the  coupling  between  the  P-SV  and  SH  waves.  Velocity 
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models  with  different  volume  heterogeneities  and  surfaee  fluetuations  are  tested  and 
eompared.  For  an  isotropie  explosion  souree  in  a  typieal  layered  model  with  high  veloeity 
top  layer,  only  weak  Lg  wave  (trapped  guided  wave)  ean  be  generated  from  the  S*  phase. 
There  is  no  energy  in  the  tangential  eomponent.  When  adding  volume  heterogeneities  or 
surfaee  fluetuations  in  the  model,  the  energy  transfers  from  leaky  mode  to  trapped  guided 
waves  as  well  as  from  P-SV  mode  to  the  SH  mode.  For  the  surfaee  seattering,  an 
exponential  random  surfaee  fluetuation  provide  stronger  eoupling  than  a  Gaussian 
random  surfaee  model.  Generally,  a  model  with  more  small  seale  topographie 
fluetuations  tends  to  generate  more  trapped  energy. 
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